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Some Monte Carlo Experiments in Computing 
Multiple Integrals 


1. Introduction. Despite the large volume of theoretical material on Monte 
Carlo methods which has appeared in the last five years, relatively little has been 
published on the experimental! side of this technique [1 ]. The numerical analyst 
who is confronted with production computation is frequently anxious to obtain 
the “feel” of the method by looking at the results of the technique when applied 
to simple problems. It is towards this end that we present here the results of 
computing the volumes of n-dimensional hyperspheres (m = 2, 3, ---,12) by 
Monte Carlo methods. These computations were carried out on the National 
Bureau of Standards Eastern Automatic Computer (SEAC). 

2. The Problem Solved: An n-dimensional hypersphere H, is defined by the 
inequality 


(1) xe t+ xe <1. 


Its volume V, is given by 


2 wn 
Va = 


(2) 


if Z,, designates the hypercube 
(3) (¢ = 1,2, m) 
then the fraction of the volume of EZ, which is contained in H, is 


The quantity », will be called the normalized volume. The ratio v, may be com- 
puted in Monte Carlo style by generating a pseudo-random sequence of points 
P,, P:, --- which lie in E, and finding the success ratio, i.e., the fraction of the 
points which lie in H,. In the computations carried out, no use was made of the 
symmetries of H,, nor was any attempt made to carry out any of the ideas of 
“importance sampling’’ [2]. 

3. Generation of Pseudo-Random Sequences: Two methods were employed 
to generate a pseudo-random sequence [3]. Method I is that of successive low 
order multiplication by 5'’-2—*. This can be described by 


(5) 
= (mod % = 1. 


which modified Fibonacci sequences. Starting from two 
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arbitrarily given values uo, u; < 4, we form successively 
(6) = + mod 4. 


The interest in the second method lies in its utilization of addition as opposed to 
multiplication and raises the possibility of consuming less time. Some preliminary 
statistical work on these sequences [4 ] indicated that it was better to select every 
second number of the sequence (6). With the aim of scrambling the fundamental 
sequence (6) even further, we introduced skips of length 1, 2, 3 in (6) depending 
upon whether the last hexidecimal digit of the last number selected was = 1, 2, or 
0, mod 3. In the runs which were tabulated under Method II, the starting values 
uo = O and u, = 2-” were selected. In the runs tabulated under Method II’, the 
starting values = 7, = were used. 

To generate a sequence of points P;, P2, ---, in E,, the following scheme was 
adopted. The fundamental sequence was normalized so that —1 < u; < 1. The 
computing machine then inquired : is u;? + < 1? Ifso, is + uo? + u;* < 1?, 
-++, etc. If all these questions were answered affirmatively, then P;: (u, ---, un) 
and P,¢eH,. If the question was answered negatively at some point, say 
+ue+--- +u,>1,p < n, then ---, were considered as the first p 
coordinates of a point P,; not in H,. The remaining coordinates up41, ---, %, were 
then not discarded, but the same procedure was repeated for P, beginning with 
tp+1. That is, if P,e H, then Pe: (tény1, Untp’, but if P, is not in 
H,, then P, has its first coordinate beginning with u,1. 

4. Quality of Results: According to statistical theory [5], the probable error 
(error at the 50% level) in a computation with N points is 


(7) eso = .6745 Sn 


while the error at the 95% level is given by 


(8) é95 = 1.96 Nee 


With Method I, the observed error was less than ég5 in 10 of the twelve cases and 
exceeded it by about 20% in two of the twelve cases. The quality of Method II 
is somewhat dubious. The error estimate ¢9, was exceeded in two cases out of six, 
and in the case m = 9, the discrepancy is pronounced. Convergence seems to be 
taking place towards the wrong value. Method II’ (which differs from Method II 
only in the starting values) seems to be considerably better and on a par with 
Method I. In view of the erratic nature of the results of Method II and its ap- 
parent sensitivity to starting values, we do not recommend that it be used for 
purposes of multiple integration. Method I seems to be satisfactory and the results 
are in conformity with the 1/ VN law for accuracy. This law is extremely con- 
servative and precludes many computations from being carried out economically, 
and if, indeed, the present speed of electronic calculators were increased by, say, 
a factor of 100 this would result in only one more decimal place of accuracy for a 
fixed amount of computing time. As it is, using Method I, the approximate 
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TABLE I 
Computation of Normalized Volume of n-Dimensional Hypersphere by Sampling 
t Method I 
No. of 
tal Points 
Sampled n=2 n=3 n=A4 n= 5 
32  0.8750000000  0.5000000000  0.2812500000  0.2187500000 
64 8281250000 .5000000000 .3125000000 .2187500000 
, OF 128 .7812500000 .5468750000 .3046875000 .2187500000 
ues 256 .7929687500 .5625000000 .3242187500 .1835937500 
the 512 .7851562500 5449218750 .3398437500 .2011718750 
1024 8027343750 5341796875 .3183593750 .1757812500 
2048 .7880859375 5307617188 .3144531250 .1708984375 
4096 .7817382813 .5219726563 .3022460938 .1638183594 
rhe 8192 .7860107422 5224609375 3067626953 1647949219 
1?, 16384 .7907714844 5299682617 .3120117188 1635742188 
Un) 32768 .7890014648 5272216797 .3111267090 1646423340 
say 
st p Exact 
aes Answer 0.785398 0.523599 0.308425 0.164493 
ith 
ve Error 0.003603 0.003623 0.002702 0.000149 
Predicted 
‘ror Error at 
95% Level 0.004445 0.005408 0.005001 0.004014 
No. of 
Points 
Sampled n=6 n=7 n=8 n=9 
32  0.1250000000  0.0000000000  0.0000000000  —0.0000000000 
64 .1562500000 .0312500000 .0156250000 .0156250000 
128 .1250000000 .0625000000 .0390625000 .0234375000 
256 .1054687500 .0468750000 0273437500 .0156250000 
512 .1113281250 .0449218750 .0312500000 0214843750 
and 1024 .0917968750 0361328125 0195312500 0126953125 
d II 2048 .0825195313 .0317382813 .0156250000 .0087890625 
six 4096 .0778808594 0319824219 0151367188 .0078125000 
be 8192 .0822753906 .0346679688 0144042969 0065917969 
> 16384 .0800170898 0342407227 .0136108398 0056762695 
d Il 32768 0792236328 0344543457 0142517090 0057983398 
vith 65536 0061340332 
ap- 
Exact 
ults Answer 0.080746 0.036912 0.015854 0.006442 
ily, Error 0.001522 0.002458 0.001602 0.000308 
say, Predicted 
ora Error at 
nate 95% Level 0.002950 0.002041 0.001352 0.000612 
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No. of 
‘Points 


Sampled 


32 

64 
128 
256 
512 
1024 
2048 
4096 
8192 
16384 
32768 
65536 


Exact 
Answer 


Error 


Predicted 
Error at 


n= 10 


0.0000000000 
-0000000000 
-0078125000 
-0039062500 
0078125000 
-0048828125 
-0024414063 
-0024414063 
.0021972656 
-0015869141 
0017700195 
0022277832 


0.002490 
0.000262 


95% Level 0.000540 


Error - 


Predicted 
Error at: 


95% Level 


n=2 


0.125 
46875 
65625 
-7421875 
-7695313 
7851563 
-7783203 
-7792969 
-7861328 
-7784424 
-7832642 
-7859802 
-7859192 


0.785398 


0.000521 


0.003143 


n= 11 


0.0000000000 
-0000000000 
.0078125000 
-0039062500 
.0039062500 
0019531250 
.0009765625 
.0009765625 
.0007324219 
.00067 13867 
-0007934570 
.0010223389 


0.000920 


0.000102 


0.000232 


Method IT 


5214844 
5380859 
5327148 


0.5235987 


0.0091 


0.021631 


TABLE I—Continued 


m = 12 


0.0000000000 
-0000000000 
-0000000000 
-0000000000 
-0019531250 
.0009765625 
-0004882813 
-0002441406 
.0001220703 
-0002441406 
.0003051758 
-0004425049 


0.000326 


0.000117 


0.000138 


-3183594 

-31591768 
-31494141 
31323242 


0.308425 


0.004807 


0.010001 


n= 5 


0.000 
.03125 
109375 
1171875 
-1406250 
-1503906 
1572266 
1645507 
1594238 


0.164493 


0.0051 


0.011353 


Number 
of Points 
Sampled n=3 n=4 
16 0.125 0.0625 
32 .3125 1875 
64 453125 .28125 
3 128 .5000000 .3046875 
256 .3125000 
512 .3203425 
1024 
2048 
4096 
8192 
16384 
32768 
65536 
Exact 
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Number 
of Points 
Sampled 


16 

32 

64 
128 
256 
512 
1024 
2048 
4096 
8192 
16384 
32768 


Exact 
Answer 


n=6 


0.000 
03125 
.078125 
0859375 
.0742188 
.0761719 
0791016 
.0786133 
-0764160 
0740967 


0.0807455 


0.0066 
0.005900 


Number 
of Points 
Sampled 


32 

64 
128 
256 
512 
1024 
2048 
4096 
8192 
16384 
32768 


Exact 
Answer 


Error 


Predicted 
Error at 
95% Level 


n=9 


0.000 
-000 
-000 
-0078125 
0039063 
0019531 
0029297 
-0053711 
0043945 
-0039063 
0040894 
0038147 


0.0064424 


0.0026 


0.000866 


Method IT’ 


n=2 


0.8125000000 
-7500000000 
-7812500000 
-7695312500 
-7714843750 
-7822265625 
-7812500000 
-7800292969 
-7800292969 


0.785398 


0.005369 


0.008890 


TABLE I—Continued 


n=9 


0.0312500000 
-0156250000 
-0156250000 
.0078125000 
.0039062500 
-0048828125 
-0048828125 
-00537 10938 
.0072021484 
.0068359375 
-0073852539 


0.006442 


0.000943 


0.000866 


5 
Error | 
Predicted 
Error at 
95% Level | 
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amount of time required to compute », is 2nNy, where yu is the multiplication time 
of the machine. The case m = 12, N = 65,536 consumed about 14 hours on SEAC. 

5. Some General Observations on Multiple Integration: One reason for 
carrying out the experiments described is this. A pseudo-random sequence is one 
which has passed (in some vague sense) a certain number of statistical tests. It is 
not clear a priori that by combining its consecutive elements into n-tuples we 
obtain a pseudo-random sequence in m dimensions. Therefore the above experi- 
ments constitute additional tests for the randomness of the sequence and, when 
successful, increase our confidence in its use in multiple integration. 

For purposes of integration, it is not necessary to employ pseudo-random 
sequences; it suffices to employ sequences which are equidistributed [6]. Such 
sequences can be defined rigorously and the integration of sufficiently smooth 
functions can be carried out with a theoretical error of order O(1/N) [7]. 

To obtain a picture of the difference in the speed of convergence, compare the 
final two tables. In the first of these tables, the integral /o' x dx was computed 
using the sequence (5). In the second table the same integral was computed using 
the sequence 


(9) uj; = E + : [x] = fractional part of x 
for which the theoretical error estimate is O(1/N). The manner in which the 
former computation obeys the 1/ VN law, and the latter the 1/N law is striking. 

To integrate in m-dimensions, one makes use of the points P;: ((jé1], Ljé2], 
--+, [jt ]) where &; are linearly independent irrational numbers belonging to 
some algebraic number field. The error estimate is again O(1/N), but of course one 
must contend with a possible increase in the size of the constant which is implicit 
in the term O(1/N). 

Table III presents the results of an integration carried out in four dimensions 
using £ = =, & = V3, & = 5 &= V10. It would appear that these figures 
satisfy the 1/N law although the associated constant is somewhat larger than that 
of the computation in Table II. 

One final observation. The use of the word ‘‘Monte Carlo”’ or “‘sampling’’ in 
connection with quadratures camouflages the following computational fact of 
life: having agreed to invest in an N point computation, we proceed to form N 
fixed points P,, ---, Py and use the approximation 


N 
(10) = V/N & f(Pi). 


This is tantamount to the use of an N point rule with abscissas P; and weights 
equal to 1/N. Such an N point rule is incredibly inaccurate. For any fixed set of 
N points, we can, in principle at least, select the weights so that the error is 
minimized in some sense [8 ]. Or we can select both abscissas and weights so that 
this error is minimized. It might prove fruitful to look into the possibility of 
obtaining best quadrature formulas of high order and in spaces of high dimension. 
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TABLE II 
1 
1e Computation of f x dx 
0 
ne Using the Sequence u; ; 2 
Described by (5) 
ve No. of Points 
Sampled 
“a 2 0.3703920880 0.5606601718 
4 5216710794 5177669530 
5064453536 5569805153 
16 4547871881 5104076401 
ch 32 .5015472909 5110118896 
th 64 5204149952 4965953886 
128 5268724155 4990123865 
256 5070362771 4999401325 
he 512 5122595896 4998424994 
ed 1024 5147580545 4996472331 
ng 2048 5140853905 4997449819 
4096 .5083722340 4999404795 
8192 5031562483 .5000873340 
16384 .5026730074 .5000148320 
32768 5019506662 .5000224160 
65536 5028558372 .5000070665 
131072 5016144215 4999992555 
the 262144 5008443891 4999988924 
ng. 
TABLE III 
1 1 1 1 
Computation of f f f f Using the Sequence 
0 0 0 0 
ans No. of Points 
Sampled 
2 1.0556385354 
hat 4 1.0646192182 
8 1.0592766300 
16 1.0615566540 
in 32 1.0626118843 
t of 64 1.0586261203 
n N 128 1.0657314237 
256 1.0673118827 
512 1.0668402647 
1024 1.0681499894 
2048 1.0685418201 
4096 1.0685544568 
8192 1.0688021121 
ghts 16384 1.0691568057 
et of 32768 1.0691964247 
or 1s Exact Answer 1.06939761 
that P. Davis 
‘y of National Bureau of Standards 
Washington, D. C. 
P. RABINOWITZ 
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Generating and Testing Pseudo Random 
Numbers on the IBM Type 701 


1. Introduction. With the increased use of large-scale digital computing 
machines many problems too complicated to solve analytically, or even nu- 
merically, have been solved by the so-called Monte Carlo method. The Monte 
Carlo technique employs repeated sampling to evaluate integrals or to simulate 
physical problems directly and requires a large supply of random numbers dis- 
tributed in a manner resembling the true distribution of the physical quantities 
which they represent. 

A set of random numbers of any specified distribution can be, and usually is, 
obtained from a set of random numbers uniformly distributed on the interval 
0 to 1. If we are given the density function f(x) of the desired distribution and 
can integrate this density function in closed form, the most direct way to get the 
desired set of random numbers {NV} from a set of uniformly distributed random 
numbers, {m;}, is to solve for NV; in the equation 


(1) f(x) de. 


However, the number of cases where the integral can be evaluated and the result- 
ing equation solved explicitly for N; is small, especially when such solution must 
be simple enough to be repeated the large number of times required for an ex- 
tensive table. In some cases the density function f(x) may be approximated by 
some simpler function, a polynomial say, to afford an approximation of (1) which 
can easily be solved for N;. 

In some other cases there are special properties of the distribution of NV; by 
which one is able to get the new set of numbers quite directly from the original 
set. One example is the normal distribution, or error function, a common distribu- 
tion of random quantities encountered in physical problems. From the central 
limit theorem of probability we know that sums of uniformly distributed random 
variables will approach the normal distribution as the number of random variables 
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being summed becomes large. In practice, adding the uniformly distributed 
random numbers in relatively small blocks, on the order of 10, produces a satis- 
factory set of normally distributed random numbers for many physical applica- 
tions. In other cases, especially where the extreme values of the distribution are 
important, such an approximation would not be satisfactory. 

When special methods are not available we can always fall back on a sampling 
method proposed by von NEUMANN [9]. This method is quite general and can 
be applied to any distribution for which the density function is known, although 
it is inefficient in some extreme cases. 

2. Generating uniformly distributed random numbers. The problem here then 
is to find a suitable method for generating uniformly distributed random numbers, 
which we shall refer to simply as random numbers throughout the remainder 
of this paper. There are tables of random digits available, notably those of 
Tippett [8], KENDALL and Situ [3], and the RAND Corporation’s digits and 
normal deviates [11]. The loading and storage of such tables within a high speed 
computer is inefficient. Much to be preferred is a scheme for automatically 
generating a sequence of random numbers within the machine and to generate 
them as they are used rather than storing a large table of such numbers. This 
limits the possibilities to a numerical method, and several schemes have been 
proposed. Many of these are described by Votaw and Rarrerty [10], and 
Taussky and Topp [7]. 

One of these methods is the mid-square method, wherein an m digit number is 
squared to produce a 2m digit number from which some middle m digits are taken 
for the next number and the process is repeated. The main drawback to this 
method is that undetected short cycles can occur. ForsyTHE [1], exhibits 
several of these degenerate cases and indicates the probability with which 
each occurs. 

In any practical scheme for generating random numbers it is impossible to 
avoid cycling. However, cycles as such are not necessarily undesirable; for a 
sequence of 10 digit numbers that cycled with a period of 108, say, would provide 
an adequate supply of random digits before any number was repeated. 

3. The method of congruences. The method of congruences, first proposed 
by LeHMeER [5], provides just this sort of sequence in which no undetected 
cycles can occur. Lehmer suggested generating a sequence {x,} from the re- 
cursion relation 


(2) Xn = mod m. 


When k = 23 and m = 108+ 1, one can obtain 5,882,352 eight digit num- 
bers before the sequence begins to repeat. It follows at once from elementary 
number theory that as long as k is chosen relatively prime to m, the first number 
in the sequence to be repeated must be the starting number x». Thus, cycling is 
easy to detect in a machine’calculation of such a sequence. 

Since the IBM 701 operates in the binary number system and carries numbers 
of 35 bits, we choose for our recursion relation 


(3) = 23x,_, mod + 1) 
xo = 10,987,654,321. 
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The reduction modulo 2° + 1 can be accomplished by dropping all bits above 
the 35th in the product and subtracting this dropped portion from the remaining 
35 bits. This can be done in a very few operations and the calculation of random 
numbers proceeds at high speed. A total of 172,900 of the 35 bit numbers were 
calculated and stored on tape in blocks of 100 in less than 20 minutes. Each 
number was checked for random machine errors by performing the calculation 
twice and was compared with x» to check for cycling. A subsequent calculation 
verified that the cycle had a period of 1,034,040, the maximum possible for the 
modulus 2*5 + 1. 

Juncosa [2] used a congruence method to generate random numbers but 
chose moduli of the form 28. Lehmer [6] points out that the low order digits are 
subject to short periods and suggested moduli of the form 2° + 1. However, 
Juncosa’s choice of a large value of k is superior to the choice k = 23 as he thus 
avoids several successive small numbers when some x is very small or several 
successive large numbers when x is very large. 

4. Analysis of random numbers. To test the randomness of a table of decimal 
digits KENDALL and Smitx [4] used four tests which they called the frequency 
test, serial test, poker test, and gap test. Because our numbers were generated 
as binary numbers and were to be used as binary numbers, we preferred to make 
the tests on the binary digits. The three tests actually performed on the binary 
digits are analogous to the first three tests proposed by Kendall and Smith. The 
gap test was not made because it could not be programmed to run on the machine 
at the same time that the other three were run. 

The problem for which the random numbers were generated required 10 bit 
random numbers; so the analysis was carried out entirely on the 10 bit pieces 
of the 35 bit numbers, two successive 35 bit numbers making seven 10 bit num- 
bers. Only the first 112,000 35 bit numbers were tested. The numbers were ana- 
lyzed in blocks of 14,000 ten bit numbers using the following three tests. 

Test I was merely a count of the number of times each of the 1024 different 
10 bit numbers appeared in a block of 14,000. In block 19, which we have taken 
as a typical block, the minimum occurrence of any number was 5, and the 
maximum was 26. 

Test II was simply a count of ones; and in block 19 there were 69,829 ones out 
of 140,000 digits. 

Test III was the so-called poker test. A count was made of each of the eleven 
different combinations of zeros and ones in a 10 bit number. The theoretically 
expected distribution was the binomial distribution. 

In all the above tests the fit of the actual distribution to the theoretically 
expected was determined by calculating x*. The results of the tests for block 19 
are shown in Table 1. 

Table 2 summarizes the results of the three tests on each of the 28 blocks of 
14,000 ten bit numbers. Certain blocks failed to pass one or more of the tests, 
notably block 13 which failed on both tests II and III. It should be noted here, 
however, that these tests are not independent; for a set of numbers which contain 
more ones than zeros would tend to be biased in favor of the combinations con- 
taining more ones on the poker test. 
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Because one or two of the blocks of 14,000 numbers fail on certain tests does 
not mean that the entire set of 28 blocks must be rejected. In fact, out of 28 values 
of P we should expect some extreme values, and might even view the entire set 
with some suspicion if none appeared. To test whether or not the number of 
extreme values of P is more than could be expected in a set of 28 such values, 
a x’ test was made on each of the three distributions of P: The results are shown 
in Table 2 and appear satisfactory for all three tests. 

In Table 3 is tabulated the total count of ones and the distribution of poker 
hands for the entire set of 392,000 ten bit numbers. The latter test yields P = .956 
showing that the fit of the actual to the expected distribution of poker hands is 
almost too good to be truly random. 


TABLE 1 
Tests for Randomness on a Typical Block of 14,000 Ten Bit Numbers 


I. Analysis of counts of the 1024 different 10 bit numbers: 
Expected count of each number: 13.67 


x? = 989 P = .1134 
II. Total number of 1’s in the 140,000 bits: 69829 
Expected number of 1’s in the 140,000 bits: 70000 
x? = .8355 P = .36 


III. Distribution of Poker Hands in the 14,000 10 bit numbers: 
Observed Expected 


10 zeros 12 13.67 
9 zeros, 1 one 149 136.72 
8 zeros, 2 ones 607 615.23 
7 zeros, 3 ones 1662 1640.63 x? = 7.372 
6 zeros, 4 ones 2922 2871.09 
5 zeros, 5 ones 3468 3445.31 = .69 
4 zeros, 6 ones 2763 2871.09 
3 zeros, 7 ones 1633 1640.63 
2 zeros, 8 ones 627 615.23 
1 zero, 9% ones 143 136.72 
10 ones 14 13.67 


5. Conclusions. One conclusion we can draw from the results presented here 
is that the numbers as a whole are probably quite satisfactory for most purposes. 

However, if a small quantity of random numbers, say 10,000, are to be used, 
we should avoid selecting them from a biased set such as block 13 in the above 
analysis. 

If anything, the results of the tests are too close to the theoretically expected. 
Of the 28 blocks we find four values of P above the 95% level and three values 
below the 5% level on test II. This is more than twice as many as should occur in 
these regions. However, the number of blocks under consideration is small and 
the result of test II for the total of all 28 blocks is satisfactory. On the other hand, 
the distribution of the poker hands for the total of 28 blocks is sufficiently close 
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to the theoretical distribution to be considered borderline as far as randomness 
is concerned. 

We should also remark here that there are more tests that could be made on 
these numbers. Some test analogous to the gap test proposed by Kendall and 
Smith should be made. Many other types of tests could also be made on these 


TABLE 2 
Result of the Three Tests on all 28 Blocks of 14,000 Ten Bit Numbers 


III. Distribution 


I. Word Count II. Count of 1’s of Poker Hands 
Block 
No. x? Count x? P x? P 
1 957 -929 69999 .0000286 >.995 14.682 as 
2 989 173 70055 .0864 .78 11.323 45 
3 1076 121 70190 1.03 .35 18.601 .05 
4 985 .799 69735 2.006 14 13.184 .23 
5 1016 .873 70018 .0093 .93 9.455 49 
6 910 .995 69987 .0048 95 8.299 .60 
7 956 -932 69750 1.786 .095 8.847 55 
8 1104 .038 69994 -00103 .975 8.431 59 
9 1031 425 70093 .2465 .64 7.437 .69 
10 1046 .302 69980 .0114 -930 13.411 .20 
11 999 .699 69802 1.1173 .29 10.723 .39 
12 1033 409 69999 .0000286 >.995 9.638 47 
13 1031 425 70666 12.641 <.005 26.345 <.005 
14 1051 .264 70202 1.163 .28 7.511 .68 
15 1029 440 69648 3.531 .065 13.143 .22 
16 1005 .652 69943 .0926 .16 7.856 .64 
17 982 .816 70230 1.508 723 13.646 19 
18 977 -846 69947 .0801 17 13.055 .23 
19 989 7173 69829 .8355 .36 7.372 .69 
20 950 .948 70405 4.675 .029 13.419 .20 
21 966 .898 70200 1.14 .29 10.046 44 
22 1090 .071 69935 .1204 713 17.013 .075 
23 1104 .038 69931 .1357 in 2.331 .992 
24 1034 .397 70017 .0082 .92- 8.553 .59 
25 1059 .212 69917 .1963 .66 12.721 .24 
26 1094 .059 69866 5117 49 2.740 -985 
27 931 .981 70002 .000114 .993 4.571 .92 
28 1036 .382 69611 4.313 .031 14.818 .125 


Distribution of x? for the Three Tests 


Expected 
Interval for P Number Test I Test II Test III 

O- .2 5.6 5 4 6 
.2- 4 5 5 8 7 
A- 6 5.6 4 1 7 
.6-— 8 5.6 5 7 5 
.8-1.0 5.6 9 8 3 

x? 2.713 6.641 2.000 

P .60 15 72 


n 
f 
e 
ti 
B 
S 
of 
Pp. 
ho 
Ma 
Math 
Mat 
— 
Mant 
Mat 
Cc 
Ma! 
Ma 
P 


GENERATING AND TESTING PSEUDO RANDOM NUMBERS 13 


numbers. The reader is referred to Taussky and Todd [7] for a description of 


other tests that could be made as well as to the bibliography of that paper for 
further information on testing. 


TABLE 3 
Tests on Entire Set of 392,000 Ten Bit Numbers 
total count of 1’s 1,960,339 
count of 1’s 1,960,000 


x? = .11726 
P = .73 


x? = 3.791 
P= .956 


Though the testing of the random numbers described here is by no means 
exhaustive, the numbers have been used successfully on Monte Carlo type solu- 
tions of problems on both digital and analog computing equipment. 


D. L. JonNnson 
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Evaluation of an Integral Arising in Numerical 
Integration Near a Logarithmic Singularity 


This paper is concerned with the evaluation of 


in terms of values of f. If this integral is known, the evaluation of 
(2) i= f ” f(x) log x dx 


is readily accomplished. For, let 


@) I, = logy 
If f is bounded, then 


v x 
= ff) de = + 


the two forms for J converge or diverge together. Equation (4) follows from (2) 
by partial integration. Now the integral in (3) can be approximated by any one 
of several standard methods. The approximate evaluation of (2) can therefore be 
achieved, once the evaluation of (1) is completed. In this note coefficients of the 
Lagrangian type are presented to facilitate the computation of J. A rigorous 
treatment of the error is also presented. 

Let 


(5) f idx f f(u)du = (rh/Dz) fj + Rar 


where the coefficients vf? /D,, are such that the remainder R,, vanishes whenever 
f(u) is a polynomial in u of degree n, 


D, = 


~ 


and D, is independent of r, and f; are the tabular values of f(u) at the points 
u = jh, 7 = 0(1)n. Thus, we assume that f(~) may be adequately represented by 
the (” + 1) point Lagrangian polynomial. The interpolating polynomials are 
easily derived and have recently been published for » = 1(1)10 by Sauzer [1]. 
The coefficients yf” and D, are then readily evaluated in a straightforward manner. 
These are given exactly in Table 1 for m = 1(1)10, r = 1(1)(# + 1). It is under- 
stood that y{” = 0 if j > m, and the m superscript is omitted where no confusion 
will result. We note some properties of the integrating coefficients useful as a 
check on the computation. 


r=j=0 
(6) =0, r=0, jf 
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As a function of r, yf” is a polynomial in r of degree and so the (nm + 1)th differ- 
ence must vanish. Hence 

(7) (r) = 0. 


Thus, employing tables of Lagrangian interpolating coefficients [2], the present 
entries are easily extrapolated. 


Treatment of the error follows the discussion given by Mine [3, 4]. Let 
c = max (nh, rh). We have 


(8) f Gur(s) (5)ds, 


where 


Thus (0,c¢) is the open interval in which G,,(s) does not vanish identically. 
Gar(s) is said to be definite if it does not change sign in this interval, otherwise 
it is indefinite. If G,,(s) is definite, the mean value theorem applies and 


(11) Rar = (8); Ene = f Gar(s)ds, 
0<éE<e. 


The latter is obtained by putting f(x) = x***/(m + 1)! in (5). If we put 
f(x) equal to the polynomial of degree m + 1 which vanishes at x = 0, h, ---, mh 
and which takes the value h"** at x = (m + 1)h and use equations (5), (11), 
and (5) with replaced by m + 1, we see that 


(13) Ene = (th/Day 


since R,+:,- = 0. Thus, except for the case m = 10, the error coefficients are avail- 
able in Table 1. The elements for m = 10 are given in Table 2. These were obtained 
by calculating +{}, and checked using (12). They are also interpolable, recognizing 


that Eyo,,/r? is a polynomial in r of degree 10. Since lim Ey,,-/r? = 1/44 if h is 
r—0 


unity, there are 12 entries, and the requirement that the eleventh difference 
vanish furnished a further check on the entries. The values were calculated to 
7D and rounded to 5D for presentation. 

The error given by (8) may also be appraised by 


(14) | < sup] f(s) | f | Gur(s) |ds. 


Alternatively R,, can be expressed in terms of an integral on one of the lower 
derivatives f™, f*-», ---, f™ of f. The only essential advantage of definiteness 
is that it simplifies the evaluation of the integral in (14). 
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TABLE 2 
Values of Error Coefficients for n = 10 


0.00666 
0.01046 
0.01257 
0.01408 
0.01524 
0.01620 
0.01701 
0.01770 
0.01834 
0.01873 
0.02454 


Concerning the definiteness of G,,(s), general theorems are wanting, but in 
each case definiteness can be checked directly. It has been established that the 
function is definite for all reported values of r and m such that r < n.Ifr = n + 1, 
Ga,n+1(S) is definite (indefinite) if m is even (odd). In the indefinite cases, Ga, 241(s) 
has one zero in the open interval (0, + 1). Let th be the vanishing point. Then 
the remainder can be decomposed into two parts. Thus, 


(15) = f° Ga, nzi(s)ds + ft (€) Ga, n4i(s)ds 


where 0 < & < {h < & < (m+ 1)h. The values of ¢ and the two integrals in 
(15) are recorded in Table 3 for m = 1(2)9, mostly to 5d. The zeros were evaluated 
by inverse interpolation using TAYLOR’s theorem. Availability of the derivatives 
facilitated the calculation of the integrals. 

The heuristic verification that G,,(s) is definite ‘is straightforward, but the 
numerical work is considerably simplified, especially for m large, using some results 
due to Barrett [5]. Let G9 (s) denote the k-th derivative of G,-(s) with respect 
to s. He has shown that G®(s) is defined and continuous in (0, c) closed for 
0 < k < n — 1. With this restriction on k, G® (0) = G®(c) = 0. If »(k) denotes 
the number of zeros of G®(s) in (0, c) closed, then »(k) > k + 2 and G,,(s) is 
definite if equality holds. If k = 0, the condition is both necessary and sufficient. 
For the problem at hand, we first examined the function G&~"(s) for m > 5. If 
v(n — 1) = n +1, then G,,(s) is definite. In a few cases, it was found that 
v(m — 1) = n + 2, and in some of these definiteness was established by showing 
that v(m — 2) = n. In the remaining cases, v(m — 2) = » + 1, and these were 


TABLE 3 
Values of Zeros of Ga,n+1(s) and Error Coefficients for the Indefinite Cases 


th (n+1)h 


— 0.08602 0.03046 
— 0.03961 0.00554 
— 0.02661 0.00214 
—0.02001 0.00098 
— 0.01602 0.00050 
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shown to be indefinite by examining G,,(s) itself. From the computational point 
of view, examination of the higher derivatives facilitates the study as the in- 
convenience of handling large powers of s is removed. 

Kaptan [6] has given coefficients to facilitate the evaluation of 


(16) e(x)de 


g(x) = x*[A (x) logx + B(x)], 1=0,1 
in terms of values of g(x). He writes 


rah e+P 

where p and p + P are number of tabular intervals 4 from the singularity to_the 
first and last of group of P + 1 consecutive ordinates. All formulas are derived 
on the basis that A(x) and B(x) can be represented by quadratic polynomials. 
For the most part, the tables can be described as follows. Coefficients c; are given 
for p = 1,5, 10 or 15 and (7, r2) = (9 + s,p9 + s + 1), 5s = 0(1)4; also for p = 1 
and (r1, 72) = (0,1), (0,5) and (0, 6). Here P = 5. If » = 0, additional entries 
are presented for p = 1, P = 4 and (r;, 72) = (s,s + 1), s = 0(1)4. The latter 
group are five term formulas and it is assumed that the constant term in the 
representation for A(x) is missing. Thus, tabular values at the origin are never 
used. All entries are to 10D. The error term is not considered. 

The principal application of Kaplan’s table occurs when A(x) and B(x) are 
not individually known, but g(x) only is tabulated. If each is known, then tables 
in present paper are more advantageous since representation of f(x) by higher 
order polynomials is possible. Present tables are useful to integrate out to ten 
tabular intervals h from the singularity after which ordinary integration formulas 
may be employed. For integration beyond 10h, recourse to Kaplan tables may 
still be appropriate depending on the nature of f(x). 

The problem treated in this note is an example of approximate product inte- 
gration, and has been studied by BEarp [7]. Other examples have been given 
by SHELDON [8] and Luxe [9, 10]. 

L. LuKe 
Midwest Research Institute 
Kansas City, Missouri 
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On Irregular Negative Determinants 
of Exponent 9n 


It may first be stated that this article is concerned only with the irregular 
negative determinants which occur in the study of quadratic forms. (The technical 
terms used below, ¢.g., ‘determinant,’ “irregular,” ‘“‘exponent,” “duplicating,” 
“opposite,” and “genera II, IV, VIII, etc.,’’ are fully explained in Mathews, 
“Theory of Numbers” (Part I., 1892) on pages 59, 178, 182, 149, 68, and 132 
respectively [1 ].) 

Only three examples of the above determinants have as yet been given in 
mathematical literature, viz., —D = 3299 (ScHotz and Taussky) [2], 6075 
(Perin) [3], and 11907 (Gauss) [4]. The writer, in his Thesis for the London 
Ph.D. degree written in 1953, gives 55 more examples, the complete table of 58 
determinants up to the limit 150,000 being: 


—D= 3299 56403 94284 117612 
6075 58563 97200 120987 
11907 60075 99387 122715 
17739 64395 101675 128547 
23571 70227 103627 132300 
24300 70956 105219 134379 
27675 75411 107163 135675 
29403 76059 107419 140211 
33075 77571 108675 140940 
35235 81675 109907 141075 
41067 81891 110700 146043 
46899 85995 111051 148419 
47628 87075 114075 ~ 149499 
52731 87723 115425 
54675 93555 116883 


A characteristic property of the determinants of exponent 9 given in the above 
list (with the exception of —3299) is that they possess 13 classes each of which 
when duplicated produces its own opposite. These may be termed “critical 
classes.” 

For the determinant —D = 6075 they are: 


7 1 1519 

9 3 676 
36 3 169 
36 15 175 
76 9 81 
79 18 81 
81 36 91 
19 9 324 
25 5 244 
31 1 196 
49 1 124 
61 5 100 
76 29 91 


‘ich 
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In the exceptional case, —D = 3299, there are only 4 critical classes, the single 
genus it possesses having 81 classes, which tabulate into 9 periods of 9 classes 
each. All the other members of the list of 58 determinants given above have 
genera II, IV, VIII, etc. 

In an attempt to find a relation between some of these determinants it was 
found that 25 of them were of the form —D = 6075 + 5832p, viz., 6075, 11097, 
17739, 23571, etc. Another series was found to be of the form 27675 + 27000), 
giving 27675, 54675, 81675, etc., [5]. 

In these series, when the critical classes (a, b, c) are arranged so that a is a 
perfect square, it is found that each member gives the same values of a and 5b; 
e.g., in the above example, 6075, the critical classes can be shown thus: 


—D = 6075 —D = 11907 —D = 17739 

a b c c c 
4 1 1519 2977 4435 
9 3 676 1324 1972 
36 3 169 331 493 
36 15 175 337 499 
81 9 76 148 220 
81 18 79 151 223 
81 36 91 forthe nexttwo 163 235 
324 9 19 members 11907, 37 55 
324 5 25 17739 the 43 61 
324 63 31 -c-values are: 49 67 
324 99 49 67 85 
324 117 61 79 97 
324 153 91 109 127 


and for all the others the c-values increase respectively in arithmetical progression. 
It is found that the determinants (with the exception of —D = 3299) can 
be arranged in 20 series as follows: 


Series 
No. Form Members 

1 6075 + 5832p 6075, 11907, 17739, 23571, 29403, 
35235, 41067, 46899, 52731, 58563, 
64395, 70227, 76059, 81891, 87723, 
93555, 99387, 105219, 111051, 
116883, 122715, 128547, 134379, 
140211, 146043. 

2 27675 + 27000p 27675, 54675, 81675, 108675, 135675. 

3 33075 + 27000p 33075, 60075, 87075, 114075, 141075. 

4 24300 + 46656p 24300, 70956, 117612. 

5 47628 + 46656p 47628, 94284, 140940. 

6 75411 + 74088) 75411, 149499. 

7 11907 + 74088) 11907, 85995. 
(11907 also appears in Ser. 1) 

8 33075 + 74088p 33075, 107163. 


(33075 also appears in Ser. 3) 
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Series 
No. 
9 
10 
11 
12 
13 
14 
15 
16 
17 
18 
19 
20 
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Form Members 


132300 + 74088000p 132300 
97200 + 1728000p 97200 
115425 + 91125p 115425 
110700 + 216000 110700 
101675 + 343000p 101675 
77571 + 474552p 77571 
56403 + 580093704p 56403 
103627 + 170584511336p 103627 
107419 + 7066834559000p 107419 
109907 + 982107784p 109907 
120987 + 64964808) 120987 
148419 + 35937000p 148419. 


In every case the factor of p (i.e., the constant difference) is the cube of the 
product of the square roots of the a-values of the 1st, 2nd, and 5th classes. 

All the critical classes of these determinants have been tested by duplication, 
and their opposites were obtained in every case. 

As a further test with a determinant beyond our assumed limit, we may take 
—D = 187596, the 4th member of Series 5, and we have: 


—D = 187596 
a 


9 

16 
144 
144 
81 
81 
81 
1296 
1296 
1296 
1296 
1296 
1296 


the last six of which reduce to: 


277 
385 
145 
217 
157 
301 


2 
3 20845 
2 11725 
301309 | 
ig 66 1333 
Pe 9 2317 
18 2320 
36 2332 
414 277 
558 385 
18 145 
306 217 
: 450 301 
137 745 
173 565 
| 18 1296 
89 901 
31 1201 | 
149 697. | 
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Choosing one of these classes, say (301, 149, 697), we can duplicate it by com- 
pounding (301, 149, 697) with (697, —149, 301). 
This is done most simply by a method devised by the writer: 


Let the classes be (a, b, ---) and (a’, b’, ---) where a and a’ are prime to 
each other, and take p and g such that either ap and a’g, or a’p and ag differ by 
unity. This can nearly always be done mentally, but when a and a’ are not small, 
the values of p and g are more quickly found as the constituents of the penultimate 
convergent of the continued fraction representing a/a’ or a’/a. The required com- 
pound class is then given by (aa’, a’bp + ab'q, ---) or by (aa’, a’bq + ab’p, ---), 
care being taken that the signs of » and g are so chosen that the smaller of the 
two products, e.g., ag and a’p, say, shall be negative. Applying this to the case 
in hand, we get: 


(697.301, 697.149(—19) + 301(—149)44, ---), 
(209797, —1973207 —1973356, ---), 

(209797, —3946563, ---), 

(209797, +39580, 7468), 

(7468, 2240, 697), 

(697, 149, 301), 

(301, —149, 697), 


which shows that (301, 149, 697) is a critical class, and each of the twelve other 
classes when similarly tested is found to be a critical class. 


R. J. Porter 
266 Pickering Road 
Hull, England 


1. G. B. MatHews, Theory of Numbers, Part I, 1892. 

A. & O. Taussky, “Die Hauptideale der kubischen Klassenkérper imaginar- 
onan Zahlkérper: ihre rechnerische timmung und ihr Einfluss auf den Klassen- 
eae Jahrbuch tiber die Fortschritte der Mathematik, 60, 1934, p. 126, J. reine angew. Math. 
171, 1934 19-41. 

3. Vows. Atti. Acad. Pont. Nuovi Lincei, 33, 1881, p. 


354-391. 
C. F. Gauss, ss Arithmeticae, Art. 306, Pvt, 1801; in Werke I, 1863, p. 371; 
by H ASER, 1889, p. 653-654. 


5. The 11th member of the second series, i.e. — 297675, has exponent 27, with 40 critical classes. 
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Selected References on Use of High-Speed Computers 
for Scientific Computation 


The author is often asked to recommend reading to orient mathematicians in 
the impact of high-speed computers on numerical analysis. The following list was 
prepared in answer to one ‘such request, but does not pretend to be definitive. 
The author is indebted to C. B. Tompxrns for several suggestions. 

For a list of books not necessarily influenced by high-speed computers, but 
highly pertinent to their use, see G. E. ForsyTHE, “‘A numerical analyst's fifteen- 
foot shelf,”” MTAC, v. 7, 1953, p. 221-228. 
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I. BOOKS 


A. D. Boots & K. H. V. Boot, Automatic Digital Calculators, Butterworths 
Scientific Publications, London, 1953, 231 p. 

B. V. BowpEn, Faster than Thought. A Symposium on Digital Computing 
Machines, Isaac Pitman, London, 1953, 416 p. 

ENGINEERING RESEARCH ASSOCIATES, INC., High-Speed Computing Devices, 
McGraw-Hill, New York, 1950, 451 p. 

D. R. HaRTREE, Numerical Analysis, Clarendon Press, Oxford, 1952, 287 p. 

Cecit Hastincs, Jr., J. T. Haywarp, & J. P. WonG, Jr., Approximations for 
Digital Computers, University Press, Princeton, 1955, 201 p. 

A. S. HousEHOLDER, Principles of Numerical Analysis, McGraw-Hill, New 
York, 1953 (extensive bibliography), 274 p. 

NATIONAL BuREAU OF STANDARDS, Monte Carlo Method, Applied Mathematics 
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AMERICAN MATHEMATICAL SocrEty, “Proceedings of symposia in applied 
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V. BarcMaNnn, D. Montcomery, & J. von NEUMANN, “Solution of linear 
systems of high order,” Institute for Advanced Study, Princeton, 1946. 
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real symmetric matrix,” ORNL 1574, Oak Ridge National Laboratory, 1954, 
107 p. 

H. H. Goipstine & J. von NEUMANN, “Planning and coding of problems for 
an electronic computing instrument,” issued in three parts as part II, vols. 1-3, 
Institute for Advanced Study, Princeton, 1947-1948. 
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HarvARD University, ‘“‘Proceedings of a symposium on large-scale digital 
calculating machinery,’”’ Cambridge, Mass., 1948. 
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F. W. J. Otver, “The evaluation of zeros of high-degree polynomials,” R. 
Soc. London, Phil. Trans., v. 244, 1952, p. 385-415. z 

E. Streret, ‘‘Relaxationsmethoden bester Strategie zur Lésung linearer 
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elliptic type,” Am. Math. Soc., Trans., v. 76, 1954, p. 92-111. 

Grorce E. ForsyTHEe 
University of California 
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Modified Quotients of Cylinder Functions 
The name in the title of this note is applied to the function, €,(z), defined by 
the following equation, 
2C,-1 (z) 
(1) G,(z) Ce) 


where C,(z) is a cylinder function [1] which satisfies the pair of recurrence 
formulae, 


(2) G(s) = Coals) + Conls) 
2C,’ (z) = C,-1(z) Cy41(2) 
€,(z) has not repeated zeros and poles with possible exception of the origin and 


satisfies the following Riccati’s equation, 


(3) +2 2m) +5=0. 


Its derivative, 


9G, (z) — 


@) dz C,7(2) 


{C,-1(2)C41(2) — 


g 
q 
Ss 
n 

Research. 
r 
4, 

or 


28 TECHNICAL NOTES AND SHORT PAPERS 


has a close relation to the famous integral due to LoMMEL, 
(5) ff = E {C(z) — | 4 


Introducing each of the three kinds of Bessel functions, J,(z), Y,(z), H,™(z) 
and H,® (z), into the equation (1) in place of C,(z), we obtain corresponding kinds 
of modified quotients, $,(z), J),(z), (z) and $,®(z), respectively. In various 
boundary value problems of mathematical physics, we encounter quite often the 
Bessel functions in quotient forms [2]. It is obvious that the modified quotients 
defined here give a convenient approach to mathematical analysis and numerical 
estimates of these problems. Moreover there is a remarkable parallelism between 
the modified quotients and the trigonometric cotangent and tangent, as there is 
between the Bessel functions and the sine and cosine. Therefore the modified 
quotients should have the same raison d’étre in cylinder functions as the cotangent 
and tangent in trigonometric functions. 

From the above consideration it seems highly desirable to give permanent 
symbols to these modified quotients, to collect their formulae and to construct 
their tables. To this end an attempt was made [3], but further cooperation and 
criticism of interested workers are necessary. 


Morio ONOE 

Institute of Industrial Science 
University of Tokyo 
Chiba City, Japan 

1. G. N. Watson, Theory of Bessel Functions, Cambridge Univ. Press, 1922, p. 82 ff. 
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Bancrort, ‘The velocity of longitudinal waves in cylindrical bars,” Phys. Rev., 

, 1941, p. 588. : 
ons ne . Hupson, “Dispersion of elastic waves in solid circular cylinders,” Phys. Rev., 63, 
» p. 46. 
H. Sunt & L. R. WALKER, “Topics in guided-wave propagation through gyromagnetic 

media,” Bell System Tech. J., 33, 1954, p. 579. 

3. M. ONok, “Formulae and tables, modified quotients of cylinder functions,” Report of the 
Institute of Industrial Science, University of Tokyo, No. 32, 1955. 


Flip-flop as Generator of Random Binary Digits 


The aim of the present note is to show that a well known electronic element 
of digital computers, the flip-flop, may be used for generating a series of random 
binary digits with equal probabilities. 

Let us consider a flip-flop as shown on fig. 1 and let A and B denote two 
possible stable states of the flip-flop. If we switch on the contact S, the flip-flop 
will be randomly set in one of its states A or B. We may obtain by the aid of the 
flip-flop a sequence of 2k random elements X;, X2, ---, Xx, (abbreviated {X2}), 
where 


ee A, if j-th switching on the contact 5S, set flip-flop in state A 
i | B, if j-th switching on the contact S, set flip-flop in state B 


and 1 < j < 2k. 


1 
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Rp=10 
R, =160 
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Fic, 1. 


In this way we may obtain a finite random series of A and B which are sta- 
tistically independent. One series produced by the aid of a flip-flop is given below: 


AABAABBABBABBBABBAABABABBABAABABABBAA 
BABABBBBBBBBBABBBABAABBB. 


Let { Y;} be the sequence of k pairs of elements of {Xa} such that Y; = X21, Xa, 
where 1 <i < k. Omitting in { Y,} all elements of the form AA and BB we 
obtain a third sequence whose elements are the pairs AB and BA only, denoted 
in the following by 0 and 1 respectively. 

Let ;(A) and »;(B) denote probabilities that j-th switching on of contact S 
set flip-flop in state A or B respectively and suppose that »;(A) and ;(B) are 
asymmetric, say ~;(A) > »;(B). Supposing that the flip-flop does not change its 
properties during two successive switchings, we may write 


(1) = poi(A) 

(2) = p2:(B). 

From 1 and 2 we have 

(3) b2i-1(A) = pais (B) - 
Because 


(4) = 


Rp } 
y (---- 
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Ve=-50v 
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and 

(5) = p;(1), 
therefore 

(6) b:(0) = 


where ;(0) and »,(1) are probabilities of zeros and ones in the i-th place of the 
sequence { Y;}. 

The procedure above described may be used for production of binary random 
numbers by automatic digital computers. In this case the manual switch S must 
be replaced by an electronic switch, of course. 


Z. PAWLAK 
Warszawa, Poland 


Numerical Solution of the Schroedinger Equation for Central Fields 


A fast program has been written for the ILLIAC to integrate the radial 
SCHROEDINGER equation 


(1) ul (r) — + + 2V(r) — Eu =0 
boundary conditions: u(0)=0 bounded 


for any well-behaved potential V(r). More generally the program can integrate 
any linear second-order differential equation which can be put in the form 


(2) — g(r)u(r) = g(r) 


with g(r) vanishing at zero and infinity, r°g(r) bounded at zero, and g(r) bounded 
for large r [1]. 

A distinctive feature of the program is the use of the Noumerov [2], [3] 
method for the integration. This is faster than techniques (such as the RUNGE- 
Kutta) which depend on making an estimate of one or more forward points and 
improving this by an iteration scheme. Here there is no iteration, but the error 
in the “estimate’’ of each forward point is of eighth order in the step-size h, so 
that the truncation error may still be kept small. The essence of the method is 
the elimination of all odd powers of h from the TAyLor expansion about any 
point by working with three points instead of two, followed by a change of 
dependent variable which removes the 4‘-term. The calculation of a forward point 
to order h® thus requires the value of the dependent variable at the six pre- 
ceding points. 

For an equation in the form (2), the required new dependent variable is 


(3) y = u — (W?/12)(gu + 
The prescription for calculating forward points is 


+ In 


1 
(4) = — + 1 — = 240 5*y,, 
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where the sixth central difference in the last term is given, to order h*, by 


This prescription can be used as soon as y is known at the first six points. 

To get this far, i.e., to start the integration, an iteration method is used: The 
difference equation (4) without the last term is solved by successive approxima- 
tions at the first six points, guessing a value of yp and improving this guess by 
applying the condition (seventh order approximation) 


(6) = 


For the inhomogeneous equation (g ~ 0), where the normalization is -not 
arbitrary, the starting value u, = u(h) is found by another iteration scheme: 
With an estimated u,, the equation is integrated out to some large r, where the 
inhomogeneity g has become negligible. Application of the condition that u be 
small and monotonically decreasing leads to an improved estimate of u;. This 
scheme, which seldom requires more than three iterations, is practical because of 
the speed of integration. (ILLIAC takes less than 40 milliseconds per h-step.) 


In its present form, the program can— 


—find eigenvalues E (permitted negative energies) 

—print out values of the wave functions u = ug at desired values of r, 
and find the nodes and extrema 

_ —evaluate matrix elements fuguz; (with ILLIAC’s electrostatic memory 

of 2° words, at 40 bits per word, there is enough space to integrate 
several wave-functions simultaneously) 

—-solve self-consistent field problems (Schroedinger equation and Porsson 
equation solved simultaneously). 


A “guide”’ to the program has been prepared, containing a detailed description 
of the mathematical method and instructions for use of the existing tapes, and 
is available in mimeographed form at the University of Illinois Digital Com- 
puter Laboratory. 

The authors are indebted to Professors JoHn Biatt and J. N. Snyper for 
helpful suggestions. 


R. A. RUBENSTEIN 
Marjorie Huse 
STEFAN MACHLUP 


of Illinois 
llinois 
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Note on Using the Reciprocal Function for a Linear Inverse Interpolation 


Let it be required to find the argument, x, by linear interpolation from a table 


of values of y(x). If there is also a table of values of the reciprocal function ei 


v(x)’ 
then a better linear interpolation may result by using the reciprocal function. 


Frequently, in applications such as in triangle or compound interest problems, 


y(x) is given as equal to a quotient B° Hence the solution of aay | for x by 


linear interpolation can be made with equal facility. For example, if tan x = is 
given for determination of x, we might just as easily determine x from the relation 

cot x = Zz? using the cot x table. It is the purpose of this note to determine the 
conditions for which interpolation will be improved by using the reciprocal 
function. 

Suppose that y(x) and z(x) are two monotonic functions such that yz = 1. 
Since the slopes of y and z have the relation, z’y + ’z = 0, and hence are always 
opposite in sign, we choose y as the increasing function. 


Now the value, x., given by linear interpolation for the argument, x, in a 
table of values of y(x) is given by 


7 
Be 


< Be < 


Xe = xX + (xe — x1) where 


(We assume that the function y(x) is positive. Roundoff errors are not con- 
sidered.) The error, Z,, due to the interpolation is given by E, = x. — x. See [1] 
for the maximum error inherent in determining a function by linear interpolation 
from tables. It is given by 


(x2 — x1) 


error| < 
lerror| < 


M, 


where M denotes the maximum absolute value of f’’(x) in the interval (x, x2). 
However, the problem presented in this paper is one of inverse interpolation. 
Hence this error formula can be applied only after the function is inverted. For 
compound interest functions, see [2]. 

The corresponding interpolated value, x,, and the error, E;, may be deter- 
mined for the reciprocal function, z. (See [3]. It should be noted that this paper 
is not considering reciprocal differences. We consider only the differences of the 
reciprocal function.) The condition that the interpolation is improved by using 
the reciprocal function is given by |E,| <|E,|. We note that E, < 0 (>0) if 
y” > 0 (<0) and E, > 0 (<0) if 2” > 0 (<0) and thus we have the follow- 
ing cases: 


| 
| 
= - 
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Case I. y” < 0 and 2” > 0, then EZ, and E, are both positive. In this case the 
result is poorer since always 


— 2 


(1) 


(x2 > —x+ (x2 — x1) = 


32 — 21 
The difference in the errors is given by 


(y — 91) (2 — y) — 
— 


Case II. y"’ > O and 2” < 0, then E, and E;, are both negative. Inequality (1) 
remains true, but now | £;| < |£,| and hence the interpolation is always improved. 

Case III. y’’ < 0 and 2” < 0. This case is impossible. The reason for this is 
as follows: From yz = 1, one obtains z’y + y’z = 0 and 2”y + 22’y’ + zy” = 0. 
Hence 2” + 2*y’’ = 2(z’)*y, and since y > 0, y” and 2” cannot both be negative. 
In fact, if y” (or 2’) is negative, then 2” (or y’’) is positive. 

Case IV. y’’ > 0 and 2” > 0, then E, is negative and E; is positive. The 
inequality is E, <|£,|, giving the condition, 


(2) E, E,= 


(3) Ey (x2 — x1) — — x) <0 
or 
(4) E, —|E,| = (x2 — x1) (1 + *) — 2(x — <0. 
Some examples follow: 
I Case Il 
2 = cosx 2’ <0 Result is always improved. 
z= cscx Result is always poorer. 
y = Vx (x > 0) y” <0 Case! 
z’’>0 Result is always poorer. 
= (x > 0) y”>0 Case IV 
E, -|E.|=0if x = 
>0 Result is improved if 
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fy = (1 +4)" forn y’>0O Case lV 
Ey —|E_| = if = 
V >0 Result is improved if 
| 
y = Sa(t) for n Case IV 
" ae Result is same as for (1 + 4)*. 
y = for n y” <0 Case! 
be C=, (n > 0) 2’ >0 Result is always poorer. 


’ 
va tan for (0 <=) Case IV 
z = coté 2” >0. 

The error in this example is given by 


tan 6 — tan 6, 


In this case | Z,| increases as @ increases from 0 tos for constant values of 6. — 6; 


and @ — 6. That this is true is seen by the relation 


tan@—tan@, sin (6 — 6;) cos 
tan @. — tan 6; sin (02 — 4;) cos 6 


Now 


Hence E; <|E£,| only if @ > 7, and in this case the interpolation is improved; 


but if @ < ; the interpolation is poorer. 


Michigan State Universi 
Michigan 


1. jj B. ScaRBorouGH, Numerical Mathematical Analysis, Johns Maphine Press, 1930, p. 107. 
2. HuGH E. Stetson, ‘‘The accuracy of linear interpolation in tables of the mathematics of 
finance,” MTAC, v. 3, 1949, p. 408-412. 
ee ied M. Mitne-Tuomson, The Calculus of Finite Differences, Macmillan and Co., 1933, 
p. 


4. W. E. Ming, Numerical Calculus, Princeton Univ. Press, 1949, p. 219-236. 


Hues E. 
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1[A, L].—A. D. Boorn, “A note on approximating polynomials for trigonometric 
functions,” MTAC, v. 9, 1955, p. 21-23. 
J,(4x) is given to 11D for m = 0(.5)12.5. 


2[B, C].—Joun von NeuMANN & Bryant TUCKERMAN, “Continued fraction 
expansion of 2#,”” MTAC, v. 9, 1955, p. 23-24. 


Let S,(x) be the sum of the first » partial quotients of the continued fraction 
for x. A table of 2 log n/log 2 and S,(2*) is given for » = 100(100)2000; the first 
expression is given to 1D; the second is given exactly. 


3[D].—Brircer Jansson, Numerisk fourieranalys fér hand-och hélkortsberékning. 
44 leaves of blue line prints from transparencies, 29.7 cm., deposited in the 
UMT Fre. 


This contains tables of sin (mm 322) and cos (on 3&2) for v = 0(1)25 and 


‘ 360 360 
m = 1(1)25 and of sin (om 382 ) and cos ( 380) for » = 0(1)100, m = 1(1)50, 
all to 4D. The tables are prepared to assist in Fourier analysis using desk calcu- 
lators or punched card machinery. A remarkably detailed description of pro- 
cedures for Fourier analysis included ; this is based mainly on an assumption that 
punched card machinery similar to International Business Machines is avail- 
able. The discussion leads up to a careful time estimate for the calculation. 


The author announces but does not include tables of sin (en 200 


cos (an 0°) for v = 0(1)300, » = 1(1)160. Tables designed for the same 


application have been prepared also by L. W. PoLLak and by Owen Mock; these 
are reviewed in MTAC, v. 5, 1951, p. 19-21 and p. 149, and in MTAC, v. 9, 
1955, p. 72 and p. 196. No earlier tables are known to contain as detailed in- 
structions for machine analysis. 


Cc. B. T. 


4[F ].—Anprew S. Anema, A table of primitive Pythagorean triangles with their 
generators, having identical perimeter totals. 13 typewritten pages deposited in 
the UMT Frnez. 

The main table is a set of 107 primitive Pythagorean triangles with identical 
perimeters. The table contains two primitive Pythagorean triangles with perimeter 
1716 (sides 364, 627, and 725; 748, 195, 772), three primitive Pythagorean 
triangles with perimeter 14280, four primitive Pythagorean triangles with 
perimeter 3 17460, and 107 primitive Pythagorean triangles with perimeter 
60850 05270 54420. 
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These were computed without aid of machines. The author, an old hand at 
lightning calculation, includes several comments about the method of calculation. 

One example from each of the sets with 3, 4, and 107 entries follows: 

From set of three: 


P = 14280 = 2-2?-3-5-7-17 
3-5-7 = 105, 2?-17 = 68 (generator a), 37 is generator b; 
A = 2a) = 5032, B = a? — B = 3255, C = a? + B = 5993. 


From set of four: 


P = 317460. 
A = 1 53868, B = 9435, C = 1 54157. 


From set of one hundred seven: 


P = 60850 05270 54420, 
A = 2763 24098 54668, 
B = 28977 68095 17195, 
C = 29109 13076 82557. 
&. 


5[F].—N. G. W. H. BEEcEr, Table of the Least Factor of the Numbers that are 
not Divisible by 2, 3, 5 of the Eleventh Million. 429 tables, 34 pages of form, 
334 cm., deposited in the UMT Fre. 


In 1949 D. JaARDEN and A. Katz published, Page 477 to D. N. Lehmer’s 
“Factor Table for the First Ten Millions” [1]. It contains the table of the least 
factor of all the numbers not divisible by 2, 3, 5, and 7, in the interval 10 017 000- 
10 038 000. I have checked this extension against the part of KuLtxK’s table, 
“‘Magnus Canon Divisorum”’ [2], and against the part of a table constructed by 
L. PoLett1 by means of his printed fasciles, Neocribrum. No discrepancy was 
found. Therefore my construction of the table of the least factors of the eleventh 
million was commenced at cycle 335, that is to say at 10 030 021. 

In the printed fasciles (conceived by L. PoLETTI), mentioned above, the least 
factors 7, 11, 13 are printed. I used ‘‘stencils” to insert the least factors 17, 19, 
23, ---, 359. All insertions in the fasciles were made by rubber stamps. Therefore 
I calculated by means of a calculating machine and using J. GLAISHER’s table of 
the differences of consecutive primes [3]: the products 367q, 373q, ---, 3313q 
(gq are suitable primes) in the eleventh million. These were then inserted in the 
fasciles. The resulting table was checked: first against KuLIK’s table, second 
against the mentioned table by L. PoLetti. The discrepancies were corrected. 
This work confirmed the primality of all the numbers of Liste des nombres premiers 
du onziéme million (plus précisément de 10 006 741 @ 10 999 997), d’aprés les 
tables manuscrites de J. PH. Kuxix, L. PoLetti, et R. J. Porter. Imprimerie 
“‘Werto,”’ Amsterdam, 1951. 


N. G. W. H. BEEGER 
Amsterdam 
Holland 


at 
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* D. ARDEN & A. Katz, Additional page (477) to D. NV. Lehmer’s Factor Table, Riveon Lemat., 
v. 


Pu. KuLix, Magnus Canon Divisorum. [See UMT 48, MTAC, v. 2, _ 139-140. 


6[F, Z].—Epcar Karst, Tables for converting 4 digit decimal fractions to periodic 
octal fractions. 4 Tables, 13 handwritten sheets, 10.8 cm., deposited in the 
UMT Fre. 


The author lists the complete period of 4, 20, 100, and 500 octal digits for the 
decimal fractions .2, .04, .008, and .0016. For any other fraction, the period is a 
cyclic permutation of one of these periods and the proper starting place is tabu- 
lated. If the fraction is not a multiple of 5, determine the first digit or two in 
the usual way (multiply by 8) and the table may now be consulted for the periodic 
digits. The tables are arranged for the convenience of the table maker, rather 
than the user. A complete description is: 


Table V—i 
5-*k Indgsk m 8*5~-‘ (mod 8) 


where k = 1(1)5‘ (k, 5) = 1m = 1(1)4-5*" i = 4(1)1 and indices are with re- 
spect to 5‘. 

J. L. SELFripGE 
University of California 
Los Angeles, California 


7(1).—Francis L. Miksa, Stirling numbers of the first kind. 27 leaves reproduced 
from typewritten manuscript on deposit in the UMT Fine. 


After completing a table [4] of Stirling numbers of the second kind, ,S,, for 
n = 1(1)(50) the writer decided to prepare a companion table of Stirling numbers 
of the first kind, S(r, 2), computed for the same range. 

Both kinds of Stirling numbers arise in formulating relationships between 
algebraic factorials and powers. For the Stirling I numbers we need only to express 
the factorials by means of power series. 

It is customary to define the coefficients in these power series as the Stirling 
numbers of the first kind, and this of course includes the algebraic sign. For our 
purpose of tabulation, however, it is more convenient to consider only the absolute 
value of the coefficients. 

The table was computed by means of the recurrence relation 


S(r,m +1) = n-S(r — 1, + S(r, 2) 
and was checked by means of 


S(r,n) = nil. 


Comparison with the table of Exact Values of the First 200 Factorials [1] showed 
that the agreement was complete. 

The most extensive earlier table known to the writer is that published during 
1900 by J. W. L. GLatsHER [2], where they are listed up to m = 23 in our notation. 
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A more detailed treatment of applications of the Stirling I numbers appears 
in C. Jorpan [3]. 


613 ing St. 


H. S. Usier, Exact value of the first 200 factorials, New Haven, Conn., MTAC, v. 1, RMT 
158, 1943-45, p. 312. 
2. J. W. GualtsHER, ‘‘Congruences relating to the sums of products of the first » numbers 
om to other sums of products,” Quart. Jn. Math., v. 31, 1900, p. 1-35. 
. C. Jorpan, Calculus of Finite Differences, Chelsea Pub. ‘o., New York, 1947, p. 142-168. 
4. Francis L. Mrxsa, Stirling Nu s of the Second Kind. Deposited in M TAC. MT Fire. 
See Review 85, MTAC, v. 9, 1955, p. 198. 


Francis L. Mrxsa 


8[1].—Kvo-cuu Ho, “Double interpolation formulae and partial derivatives in 
terms of finite differences,” MTAC, v. 9, 1955, p. 52-62. 


Expressions for partial derivatives of h®k*-*(d*f/ax®dy**) where a = 1(1)4, 
8 = 0(1)a@ are given in terms of the differences of f at intervals # in x, k in y. 
The expressions given include the fourth order terms. These are given in various 
forms suitable for use at the edges and corners of a table as well as in the center. 
The corresponding Lagrangian expressions are given for a = 1(1)3. 


o[1}.—H. E. Sauzer, ‘‘Osculatory quadrature formulas,” J. Math. Physics, 34, 
1955, p. 103-112. 


Quadrature formulas such as 


I(r, s) = f(x)dx =hD AS + Refer (6) 


can be found by integration of osculatory interpolation formulas. The exact values 
of the A,, B;, R, expressed as rational fractions in their lowest terms, are given in 
the case of m points for m = 2(1)7. For each m the appropriate coefficients are 
given for all possible sub-intervals: thus for » = 4, with f and f’ given at points 
x; = Xo + th, i = — 1, 0, 1, 2, we find expressions for the integrals I(r, s) for 
the following values of (7, s): (—1,0), (0,1), (1, 2); (—1, 1), (0, 2); (—1, 2) in 
terms of all the data. 3 


10[1].—T. N. E. Grevitte & H. VAuGHAN, “Polynomial interpolation in terms 
of symbolic operators,” Trans. Soc. Actuar., v. 6, 1954, p. 413-476. 


Table 6, p. 450-451, gives the exact values of the basic function L(x) = L(—x) 
of the following continuous operators: M*, M*, #M*, 
M®, for x = 0(.1)4 or less (usually around 3), up to the 
value of x where L(x) = 0. The basic function L(¢) for the continuous operator 
K is defined by K = /_% L(t)E~-‘dt where E denotes the displacement operator 
of the calculus of finite differences. The operator K is also called the characteristic 
operator for the interpolation formula of which L(t) is the basic function, and it is 
analogous to the graduation operator which the authors have associated with a 
discrete interpolation formula. 
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Preceding Table 6, on p. 443-449 are the following collections of formulas: 

Table 1, p. 443-445, Characteristic Operators of Certain Published Interpolation 
Formulas, giving name or originator of formula, place and date of publication, 
and characteristic operators. 

Table 2, p. 446, gives the symbolic expansions of certain operators, 
namely, combinations or separate powers of M, yw, and 6 in terms of 


Das. (6 = Et — E4, yp = 3(E) + E+), and M is a new operator defined by 


Mf(x) = fe + 

Table 3, p. 446, Traces of powers M* and mean powers »M*, k = 1(1)8, in 
terms of even powers of 6, where the author defines the trace of a continuous oper- 
ator K with basic function L(#) as the discrete operator t(K) = >°.. L(n)E-. 

Table 4, p. 447, Powers and mean powers of M, up to the 10th, in terms of 
special operators for both endpoint and midpoint formulas. 

Table 5, p. 448-449, Various special operators for endpoint, midpoint, or 
mixed endpoint-midpoint formulas in terms of M, 6, uw, and D. 


H. E. 
Diamond Ordnance Fuze Laboratories 
Washington 25, D. C. 


11[K].—Tue RAND Corp., One Million Random Digits and 100,000 Normal 
Deviates, The Free Press, Glencoe, Illinois, 1955, 26.5 cm., xxv + 200 p. 
$10.00. The tables are also available on punched cards from The RAND 
Corporation, 1700 Main St., Santa Monica, Calif. $25.00 an order +$3.00 
per 1000 cards. . 


This book contains one million random decimal digits and 100,000 normal 
deviates, computed from the random digits to 3D. This succeeds [1 ] as the largest 
table of randon. numbers known to be available. The book also contains intro- 
ductory material describing production of the digits, tests of the digits, production 
and tests of the normal deviates, and use of the tables. 

Random numbers have long been used in various ways for the simulation of 
random or imponderable effects in various calculations. In recent years a growing 
popularity of the so-called Monte Carlo type of calculation in various shielding 
problems in physics and other problems of mathematics and applied mathematics 
has led to a greatly increased demand for random numbers; some of this recent 
work is described in [2], [3], [4]. The increasing power of modern computing 
equipment has also led to higher demands for random numbers in their more 
classical uses. The present publication was prepared in order to meet these in- 
creased demands. The largest earlier publication [1] is reported to have been 
completely inadequate in size for the needs of The RAND Corporation. 

Generally speaking, workers have been content with the generation of numbers 
conforming with frequency distribution requirements and believed to be essen- 
tially independent of the process being simulated in their method of production. 
Many systematic ways of generating such numbers have been proposed for digital 
computers; references to some of these are contained in [6], and other schemes 
may be found in [7]. The present numbers were not generated by any known 
formula, but they were rather generated by a process which (if current physical 
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theory is correct) was completely random, yielding output which could not have 
been predicted with any amount of knowledge available before the generation 
of the numbers. The digits were prepared by driving electronic counters by elec- 
tronically generated noises. (Another random method of generating binary digits 
is described earlier in this issue [11 ].) 

The random digits were subjected to numerous tests described in the intro- 
ductory material. The originally generated digits showed unacceptable biases in 
their distribution, and a new table was generated systematically by adding pairs 
of digits modulo 10 to improve (that is, to flatten) the distribution. This was 
successful. Standard tests of frequency distributions, poker hands, frequencies of 
digraphs (ordered pairs of digits), and lengths of runs of a single digit were carried 
out, presumably on the digits as they were punched on the IBM cards. 

The printed version was prepared automatically from these IBM cards, and 
it is unlikely to differ from the punched card version in more than a few digits. 

All the standard tests as reported were passed adequately by the digits. 
However, one test was improperly analyzed, and the authors have requested that 
the following correction be made in the text. It is believed that the desirability 
of some such modification was brought to the Corporation’s attention by I. J. 
Goop; it is based upon his analysis of a serial test [8]. 

Replace the paragraph beginning at the bottom of page xv by the following: 


Table 5 can be tested by a criterion originally due to Kendall and Smith and 
revised by Good. Assuming all pairs equally likely we get a normalized sum of squared 
deviations of 107.8. However, this statistic does not have a x?-distribution. On the other 
hand, it is the sum of the error variation and twice the row (or column) variation, 
where under the assumption of perfect randomness, the error variation 1s asymptotically 
distributed like x? with 81 degrees of freedom. We take the error variation as our test 
criterion. This gives a x? of 107.8 — 2(7.56) = 92.7, which is about the 0.18 level 
for 81 degrees of freedom. 


The modification above cannot be interpreted as challenging the randomness 
of the digits. 

It seemed conceivable that the standard tests reported as having been made 
in the introduction to the book would have failed to detect variable instability 
of counters (which, in any event, would be unlikely to have a bad effect on the 
randomness of the digits after the smoothing transformation made above). The 
reviewer did test a hypothesis that such instability occurred and led to excessive 
lengths of coupon collector’s sequences [9]. Accordingly the digits were trans- 
formed from decimally punched card form to binary coded decimal form suitable 
for introduction to the SWAC computing machine located at the University of 
California at Los Angeles. The coupon collector’s test described by GREENWOOD 
was run by the reviewer with the help of F. H. HoLtanper, J. L. SELFRIDGE, 
and Davin A. Pope. The longest single sequence which was required in order to 
provide at least one of each of the 10 digits had length 115 digits; the total number 
of sequences tested before the digits were exhausted was 34,248; the longest 
sequence of length 115 was not inordinately long. Based on the total number of 
sequences observed and the probabilities reported by Greenwood [9], the follow- 
ing comparison is possible between the numbers of sequences observed of various 
lengths and the numbers expected. 
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ve Length of Number Number 
on Run Observed Expected 
- 
its 10 14 12.428 
11 63 55.926 
- 12 152 143.542 
= 13 284 276.832 
Irs 14 462 446.786 
as 15 640 638.190 
of 16 813 834.267 
ed 17 1017 1020.064 
18 1264 1184.234 
nd 19 1342 1319.539 
ts. 20 1415 1422.517 
ts. oe 1480 1492.752 
iat 22 1528 1532.020 
ity 23 1515 1543.501 
J. 24 1576 1531.128 
25 1488 1499.087 
26 1431 1451.483 
nd 27 1375 1392.115 
ait 28 1299 1324.359 
™ 29 1189 1251.116 
30 1269 1174.806 
lly 31 1136 1097.397 
est 32 977 1020.446 
wal 33 972 945.150 
34 839 872.401 
35 792 802.832 
36 689 736.868 
de 37 695 674.765 
38 618 616.647 
39 554 562.532 
~ 40 494 512.363 
a 41 490 466.026 
te 42 408 423.364 
be 
45 338 315.534 
46 304 285.622 
47 265 258.376 
bey 48 235 233.594 
eat 49 196 211.081 
of 50 196 190.651 
ow- 51 188 172.129 
ous 52 148 155.352 
53 139 140.166 
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Length of Number Number 
Run Observed Expected 
Vi 

54 114 126.428 
55 117 114.009 
56 104 102.788 
57 84 92.652 
58 74 83.502 
59 70 75.243 
60 67 67.793 
61 53 61.072 
62 62 55.012 
63 53 49.549 
64 42 44.624 
65 47 40.186 
66 37 36.186 
67 22 32.583 
68 41 29.337 
69 21 26.413 
70 28 23.780 
71 20 21.408 
72 15 19.273 
73 18 17.349 
74 14 15.618 
75 11 14.059 
76 and more 102 126.631 


The reviewer cannot find any significant indication in these data to support a 
hypothesis that the digits are not random. A calculation of x? based on the as- 
sumption that 34,248 samples were divided into 67 categories with 66 degrees of 
freedom gave a value of x? about 0.1350 below the expected mean value. More 
detailed results of this coupon collector’s test may be obtained by addressing 
the reviewer. 

The normal deviates were computed from the random digits according to a 
standard procedure, see [6]. In particular, if D is a five-digit number taken from 
the table of random digits and read as an integer, then the numbers listed are 
values of x solving the equation 


(D + 0.5)10-5 = dt, 


For the right hand number of this equation the National Bureau of Standards 
tables of the normal distribution [10] were used. 

The printing is satisfactory. Reproduction was by photographic offset printing 
from pages printed by an International Business Machine Model 856 Cardatype. 
Each page of the table of random digits contains fifty lines; each line contains a 
serial number and ten groups of five decimal digits. The digits are easily legible. 
Each page of the table of normal deviates contains fifty lines; each line contains 
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a serial number and ten deviates. The material was proofread only on a sample 
basis, but the probability of serious error in printing seems low. 
C. B. T. 


1. M. G. Kenpa.t & B. B. Smitn, Random Sampling Numbers, Cambridge University Press, 


. National Bureau or STtanparps, The Monte Carlo Method ( of a Symposium 
hed in 1949), NBS AMS 12, U. S. <4 Printing Office, Washington, D. C., 1951. 


3. J. H. Curtiss, “‘Sampling methods applied to differential and difference equations,” Seminar 
on Scientific Computation, nternational Business Machines Corp., New York, 1949. 

. H. Kann, » Rageestene of Monte Carlo, The RAND Corp. (to be published). 

P. Davis & P. RaBinow!1z, “Some Monte Carlo experiments in computing multiple in- 
tara MTAC (this issue), v. X, 19 


D. L. Jounson, “‘Generating and en tis 701,” 
MTAC (this issue), 1986. 


7. L. Sacco, Manuel de Crypopaphi Payot, Paris, 1951, Chap. XI, p. 74-119. 
eee Good, “The serial test or sampling numbers and other tests for randomness,” Camb. 
Phil. Soc., Proc., v. 49, 1953, p. 276-284. 


9. Robert E. GREENWOOD, ‘ ‘Coupon collector’s test for random digits,”” MTAC, v. IX, 1955, 
p. my (Note also Corrigendum, MTAC, v. IX, 1955, p. 229.) 


1 ables of Probability Functions, 2nd ed.. U. S. Gov. Printing Office, Washington, D. C., 


Pawtak, “Flip-flop as Generator of Random Binary Digits,” MTAC (this issue), 
v. x. 1956, p. 28. 


12[K].—NBS Applied Mathematics Series, No. 44, Table of Salvo Kill Prob- 
abilities for Square Targets, U. S. Gov. Printing Office, Washington, D. C., 
1954, ix + 33 p., 26 cm. Price $0.30. 


Let {x, y, 0} denote the circular normal distribution centered at (x, y) with 
standard deviation o in each direction. A salvo of N bombs is “‘centered”’ at (£, 9), 
which has the distribution {0, yo, 4}. The bombs of the salvo are independently 
distributed according to {£, 7, cz}. If a bomb hits the square target (x* < 1, 
y? S 1), there is chance Px of a kill. Assuming the bombs act independently, 
the chance Px that the salvo kills the target is computed for Px = .1, .4, .7, 1; 
yo = 0, 1, 2, 4, 7, 11, 16, 22; o4, on = 1, 2, 4, 7, 11, 16, 22; NW = 1, 5, 10, 25, 50, 
100, 150, 200. The entries are 4D with possible error of two in the last place. 
In an introduction, A. D. HESTENEs warns the user against interpolating in this 
quintuple-entry table, except in the Px direction for small N. 


J. L. Hopcgs, jr. 
University of California 
Berkeley, California 


13(K].—K. C. S. Puttar & K. V. RAMACHANDRAN, “On the distribution of the 
ratio of the ith observation in an ordered sample from a normal population 
to an independent estimate of the standard deviation,’’ Ann. Math. Stat., 
v. 25, 1955, p. 565-572. 


Tables to 2D are given for (1), the 95th percentile of distribution of g, = x,/x; 
(2), the 5th, and (3), the 95th percentiles of the distribution of u, =|x,/s| where 
x, is the largest of a sample of m observations from a normal population with 
zero mean and unit variance and s is independently distributed with » degrees 
of freedom, for (1): 2 = 1(1)8: » = 3(1)10(2)20, 24, 30, 40, 60, 120, ~; (2): 

= 1(1)10: » = 1(1)5(5)20, 24, 30, 40, 60, 120, © ; (3): m = 1(1)8: » = 5(5)20, 
24, 30, 40, 60, 120, «, respectively. 
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A table used in the derivation of the above tables is given to 8S for coefficiesits 
a;™ which are defined by 


- 


for « = 0(1)30; k = 1(1)7. 


University of Californi: 
Loe Angeles, California, 
14[K].—P. B. Patnaik, “‘A test of significance of a difference between two sample 


proportions when the proportions are very small,” Senkhyd, v. 14, 1954, 
p. 187-202. 


The paper is concerned with the following problem: “if there is a random 
sample of N, individuals from a population, of which x; have the characteristic A, 
and a random sample of WN, from a second population, of which xz have A, then 
it is desired to test whether the chance of possessing A is the same in the two 
populations.”’ Consideration is restricted to cases in which Poisson distributions 
are acceptable approximations to the binomial distributions which give the 
probabilities for various values of X; and of X2. The problem thus becomes one 
of testing, using one observation from each Poisson distribution, the hypothesis 
that the ratio of the parameters of the two Poisson distributions is a given con- 
stant. Two test procedures are proposed. Randomized tests are not considered. 
(A) It is pointed out that, if the conditional probabilities of falling in the critical 
region for fixed sum of the two observations is kept below some preassigned value, 
the size of the critical region is, for small sum of the parameters, very much below 
this preassigned value. One of the many possible methods of interpolating is pro- 
posed to keep conditional probabilities near, but not necessarily below, the pre- 
assigned number and to make the size of the critical region near this number. 
Some tables illustrate the effects on conditional and unconditional sizes of critical 
regions. (B) The idea of a “test with minimal bias” is introduced. The power 
curve must be at or below the nominal size at the parameter value for the hy- 
pothesis under test, derivatives at this point must satisfy certain conditions, 
and among tests satisfying these conditions, that one is chosen which minimizes 
the length of the interval of parameter values for which the power curve is below 
the nominal size. The requirement of minimal bias is imposed on the conditional 
regions, making the test that of a binomial variate. Critical values for a binomial 
variate are given (in Table 6) for a(nominal size) = 0.05, 0.10; (parameter of 
the binomial) = 4, 3, 3, 2, m(index at the binomial) = 6(1)25. Several typo- 
graphical errors in titles are obvious. 


W. J. Drxon 


K. J. ARNOLD 
University 


15[K ].—Benjamin EpstTEIn, “Truncated life tests in the exponential case,”’ Ann. 
Math. Stat., v. 25, 1954, p. 555-564. 


With m items on a life test, it is decided in advance that the experiment is 
terminated either at a fixed truncation time Ty or at the time X (ro, m) necessary 
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for the roth failure, whichever is smaller (stopping rule), and it is assumed that 
the distribution of lives X(X > U) is exponential with unknown mean #6. If 
X (ro, m) < To(> To) the null hypothesis @ = 4 is rejected (accepted). Two cases 
are considered, non-replacement and replacement from the same population. 
Formulae are derived for Es(r), the expected number of observations, and E,(T), 
the expected waiting time necessary to reach a decision. In the replacement case 


they are related simply by E,(T) = ~ E,(). The probability L(@) of accepting 


6 = 60 where @ is true is calculated as the probability of reaching a decision re- 
quiring up to r failures. However, in Epste1n’s formulae (2), (13), (21), the 
sign = should read >. 

Table 1 gives the values of r and x*,_.(2r)/(2r) such that the test based on 
using 0,., > C = Oox*:-2(2r)/2r as acceptance region for @ = 6 will have the 
probabilities L(@)) = 1 — a and L(@,) > 8 for the errors of the first type with 
a = .01, .05, .10 and of the second type with 6 = .01, .05, .10 for 6./0, = 3, 2, %, 
3, 4, 5, 10. (The heading of the last column should be 8 = .10 instead of 8 = .01.) 
Table 2 gives some integer values of m = [@ox*:.(2re)/2T.] to be used in trun- 
cated non-replacement procedures for a = .01, .05; 8 = .01, .05 for 6/@, = 2, 
3, 5 and o/T> = 3, 5, 10, 20. 

Excellent practical illustrations are given as solutions of problems of the 
following type: Find a truncated replacement plan for which Ty = 500 hours 
which will accept a lot with mean life @ = 10,000 hours at least 95% of the time 
and reject a lot with mean life @ = 2000 hours at least 95% of the time. Then 
a = B = .05. The values L(6), E,(T) and E,(r) are computed for both values of 8. 


E. J. GUMBEL 


Columbia University 
New York, N. Y. 


16[K].—D. B. Duncan, ‘‘Multiple range and multiple F tests,” Biometrics, 
v. 11, 1955, p. 1-42. 


Let g(, m2) = w/s where w is the range of » independent normal variables 
having the same mean and unit standard deviation, and ms? is distributed inde- 
pendently of w as chi square with m2 degrees of freedom. 

Table II labeled ‘Significant studentized ranges for a 5% level new multiple 
range test” lists the (.95)”— percentiles of the sampling distributions of ¢(p, m2) 
for p = 2(1)10(2)20, 50, 100 and mz = 1(1)20(2)30, 40, 60, 100, «©. Table III 
lists the (.99)?— percentiles for the same distributions. These tables are used in 
performing comparisons among the means of equal sized samples from p popula- 
tions. Tables of the 95th and 99th (same percentiles for all ») percentiles of 
q(p, m2) are given in Biometrika Tables [1]. 

Care should be taken not to confuse the 5% level heading in the table (what 
the author calls 95% two mean protection level) with 5% level of significance in 
an analysis of variance test since they actually refer to different probabilities. If 
an analysis of variance test is to be performed at the 5% level to test the hy- 
pothesis that the p populations have equal means then, if the p means are equal, 
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there is a 5% probability that the data will be such that the hypothesis is re- 
jected. For the tests suggested in this paper there is a 5% least upper bound to 
the probability that any one pre-specified pair of populations having equal means 
will be judged to have unequal means. Using the suggested test for a case where 
all » populations have equal means the probability that at least one pair will be 
declared significantly different is 1 — (.95)?—. 

This is an attempt to make the suggested tests more comparable in overall 
level of significance to the case in analysis of variance where single degrees of 
freedom are used, and repeated tests at the 5% level are used. 

Tables I and IV describe an example. 

Table V classifies several discussed test procedures according to the measure 
of variability used among the means (range or variance), whether having other 
observed means between two specified means changes the difference labeled 
significant, and whether these differences depend on a fixed level a or a level which 
changes with the number of means between the two specified ones. 

Table VI lists for comparison some ‘“‘significant ranges for 5% level tests’’ of 
four discussed tests. Table VII compares the power of two classification of tests 
(p = 20, 5% level) to recognize a difference in the means of two pre-chosen 
populations. 

Due to the fact that ‘5% level’’ means different things for the different tests 
Tables VI and VII are likely to give a somewhat slanted comparison of the various 
discussed tests. 

FRANK MASSEY 


Eugene, 
1. E.S. Pearson & H. O. Hartley, Biometrika Tables for Statisticians, v. 1, Cambridge, 1954. 


17[K].—E. C. Fretter & H. O. Hartiey, “Sampling with control variables,” 
Biometrika, v. 41, 1954, p. 494-501. 


Let x and y be jointly distributed, and let f(x, 6), g(y, @) be the respective 
marginal distributions. Information about the control variable x is available, 
e.g., f(x, 0) may be completely known or @ may be known. The authors are con- 
cerned with using this information in the estimation of ¢@ for the case where 
are observations from an N(u,1) population, x = — X)?, 
y = Xmax — Xmin; the distribution of y is assumed unknown and that of x known. 
In the determination of the variance of the proposed estimate, the expectation of 
the reciprocal of a binomial variate, excluding the zero class, is needed. 


Values of (1) E'(np, n) = 5 +(*) for m = np = 1(1)10, n = 25, 
r=] 
50, © are computed directly. For m = 12(2)20, m = 25, 50, 100, 
1 1 2! 
m+p* (m+pym+2p) (m+ 


is used; for m = 25(5)45, m = 50, 100 a control is used with (2). E’(m, «) may 
be calculated from E’(m, ~) = e~™{Ei(m) — log. m — y} using available tables 


(2) E’ (m,n) = 


J 

bi 
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[1, For large mp, the approximation E’ ~ —— (mn — 


is given. Values of (1) for m = 1(1)10 are given to 4D, all others to 5D. 
INGRAM OLKIN 


Mi State University, East Lansing, Mich. 
win niversity of Chicago, Chicago, IIl. 


1. BAAS Math. Tables, v. 1, London, 1931. 
2. E. JAHNKE & F. Empe, Tables of Higher Functions, Sth ed. Leipzig: 1952 


18[K].—P. N. SomErvILLE, “Some problems of optimum sampling,” Biometrika, 
v. 41, 1954, p. 420-429. 


Let 
FO) = = (24) exp — tr yE-y’, 


where y = (y1, = (04s), = 1, = and let F(x) = Se f(y) 
where R: — © <9; <x, 1 =1,---,k. In Table 3.1 the author computes 
F,(x), x = 0(.1)2(.5)3 to 5D for k = 1, 2; 4D for k = 3,4; 3D for k = 5. F,(x) 
is available in a number of tables. F;(x) was computed using Table VIII of part 
II of the Pearson Tables [1], and F,(x). For F,(x), k = 3, 4, 5, an expansion 
of f(y) in terms of Hermite polynomials was used [2]. 


INGRAM OLKIN 
an Ui State University, East Lansing, Mich., 
and University of Chicago, Chicago, III. 


. K. PEARSON, Tables for Statisticians and Biometricians, London, 1931, part II. 
. W. F. Krsste, “An extension of a theorem of MEHLER’s on Hermite polynomials.” Camb. 
Phil: Soc., Proc., v. 41, 1945, p. 12-15. See p. 15. 


19[K].—H. F. Donce, “Skip-lot sampling plan,” Industrial Quality Control, 

v. XI, no. 5, 1955, p. 3-5. 

For derivations, formulae, and basic procedure, see the author’s earlier 
paper [1]. The procedure developed there is as follows: (1) inspect 100% of 
units consecutively as produced until i successive good units are found. (each unit 
is good or bad), (2) thereupon, inspect only an unbiased proportion f, (3) if in 
the proportion f a bad unit is found, revert to 100% inspection until 7 successive 
good units are found, etc., (4) replace bad by good units. For given f and i, the 
average outgoing proportion defective (AOQ) as a function of the incoming pro- 
portion defective (p), and (consequently) the maximum AOQ (AOQL) are 
determined. 

In the current paper, the series of units are replaced by series of lots, inspection 
of a unit by analysis of a sample, bad unit by bad lot (lot whose sample fails to 
meet a standard). Corresponding changes of meanings of p, AOQ, and AOQL to 
apply to lots rather than units. 

Figure 1 gives the AOQL = .01(.01).06, .08, .10 as a function of ¢ = 1(1)50, 
(2)100 and f = .01(.001).05(.002).10(.01).50. Figure 1 is essentially an abridge- 
ment of Figure 3 of his 1943 paper cited above (p. 272) though the AOQL curves 
are now fully shown for low 7. 


HAROLD FREEMAN 
Massachusetts Institute of Technology 
Cambridge, Mass. 


PA. Tg Donce, ‘“‘A sampling plan for continuous production,” Ann. Math. Stat., v. 14, 1943, 
Pp. 
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20[K ].—J. Arrcutson & J. A. C. Brown, “‘An estimation problem in quantitative 
assay,” Biometrika, v. 41, 1954, p. 338-343. 


A popular model for the quantitative response u to a stimulus of concentration 
xis:u = HP(a + Bx) + ¢, where H is the maximum expected response ; P(a + Bx) 
is the normal cumulative based on a linear function of x; and ¢ is N(0, o”), inde- 
pendent of x. The authors consider a model of the form: 

Inu +e=—=K+inP+e, 
t.e., Var (u) is proportional to E(u) }. Iterative equations, familiar to probit 
users, are set up to estimate a, 8, and K. 

A table of working probits, to facilitate the solution of these equations, is 
given for initial guessed probits, Y = ao + box = 1.0(0.1)9.0. This table contains 
values of the minimum working probit to 4D; of the auxiliary variable, P/Z 
to 4D or 5S; and of the weighting factor, Z*?/P? to 4D. The working probit for 


this model is y = + 
R. L. ANDERSON 


North Carolina State College 
Raleigh, N. C. 


21[K].—A. R. JoncKHEERE, “‘A distribution-free K-sample test against ordered 
alternatives,” Biometrika, v. 41, 1954, p. 133-145. 


Let X;,7 = 1, 2, ---, k, be independent random variables with the continuous 


cumulative probability functions F;(x), and let X;;, 7 = 1, 2, ---, m:, be samples 
of size m; of the X;. To test the hypothesis that X;, ---, X; have the same dis- 
tribution (Ho): Fi(x) = F2(x) = --- = F,(x) against the alternative that these 


random variables are stochastically increasing (A): Fi(x) > F2(x) >--->F;(x) 
(on p. 134 the sense of these inequalities is inverted, as it is in the footnote on p. 
135), the author considers the following statistic: Let p;; = number of pairs X;,, 
X;, such that X;, < among all possible mn; pairs, for given 7 1, k—-1; 
j =it+i1,---,k. The test statistic is S=2 mam, 
i<j i=1 
identical with a statistic considered by M. G. KENDALL [1] and closely related 
to the rank correlation coefficient. 

In his study of the probability distribution of S under (Ho) the author obtains 
the first four cumulants, studies their extreme values, and obtains asymptotic 
distributions under several assumptions. He finds in particular that if at least 
two of the sample sizes m; tend to infinity, then the limit distribution of S is 
normal, and he proposes a somewhat better approximation by STUDENT’s dis- 
tribution. Finally, the exact distribution for small samples is given. 

Table 3 gives Pr {S > So} to 3 or 4S for k samples, each of size m, covering 
the range: k = 3and m = 2,3,4,5;k =4andm = 2,3,4;k = Sandm = 2,3; 
k = 6 and m = 2; So = 0(2)--- until Pr {S > So} becomes negligible. 

Z. W. BIRNBAUM 


1. M. G. KENDALL, Rank Correlation Methods, London, 1948. 
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22[K].—R. A. Brap.ey, ‘Rank analysis of incomplete block designs II. Addi- 
tional tables for the method of paired comparisons,’ Biometrika, v. 41, 1954, 
p. 502-537. 


The present article extends the 1952 tables of BRADLEY and Terry [1] to 
include 7 and 8 repetitions for 4 treatments (¢ = 4; = 7,8) and 1 to 5 repetitions 
for 5 treatments (¢ = 5; m = 1, 2, 3, 4, 5). The problem considered is to estimate 
true treatment ratings ---,#:; 2a; = 1) and test the hypothesis all 
x; = 1/t. The tables give, for all combinations of total rank sums, estimates of the 
x; to 2D, the value of the testing statistic, B, to 3D, and the significance prob- 

t 


ability for B, to 4D. Only two treatments are tested at a time, giving 2 pairs 
per repetition. The favored treatment in each pair is given a rank of 2 and the 
other a rank of 1; hence, the total rank sum for any treatment can fall between 
n(t — 1) and 2n(t — 1). 

The author includes (i) a study of the usefulness of an approximate statistic, 
(ii) a discussion of various uses of B,, and (iii) two errata for the original tables. 


R. L. ANDERSON 
North State College 


R.A. Braptey & M. E. Terry, “Rank analysis of incom block designs. I. The method 
Biometrika, v. 39, 1952, p. 324-345 [MTAC, v. 8, 1954, p. 


23[K].—Epwarp Wa ter, “Uber die Ausnutzung der Irrtumswahrschein- 
lichkeit,”’ Mitteilungsblatt fiir Math. Stat., v. 6, 1954, p. 170-179. 


Tests of significance based on statistics which can assume only a finite number 
of different values usually do not attain their asserted level of significance. It is 
always possible to modify the test (by introducing a randomized test) in order to 
attain the asserted level of significance. The author is interested, however, in a 
modification based on the observations. He considers the sign test for testing the 
null hypothesis that the distribution function is symmetric about .zero. The 
modified sign test is defined as follows. Let y equal the number of negative ob- 
servations in the sample of m; let y’ be the rank of the negative observation with 
largest absolute value among the absolute values of all the observations; let y”’ 
be the rank of the negative observation with second largest absolute value among 
the absolute values of all the observations; etc. The null hypothesis is to be re- 
jected if y < yo, or y = yo and y’ < yo, or y = Yo, y’ = yo’ and y” < yo”, etc. 
In table 1 values of yo, yo’, yo’, --- are given for m = 5(1)25, a = .01, .05, and 
for one- and two-sided tests. The two-sided test is the same as the one-sided test 
except that one considers either the number of positive observations or the 
number of negative observations, according to which is the smaller. 

In table 2, empirical results are tabulated to compare the power of the sign 
test, modified sign test, and the #-test for a = .01, .05 in the case where the dis- 
tribution is normal with variance 1 and mean yu = .6, 1, 1.4. 

In table 3, the effective levels of significance of the sign test and modified sign 
test are given for m = 5(1)30. 


Cyrus DERMAN 


Columbia University 
New York, N. Y. 
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24[K ].—Sxozo Sumapa, “Power of R-charts,’’ Reports of Statistical Application 
Research, Union of Japanese Scientists and Engineers, v. 3, 1954, p. 70-74. 


Let x1, X2, x3, x4 be independent and normally distributed random variables 
with means m,, m2, m3, m, and common variance o*. The range R is defined as 
the difference between the largest and smallest value among the four successive 
observations. 

If m, = mz = m; = my, then the range provides an estimate of o. On the 
other hand, if not all the m,’s are equal then the range may include the effect of 
the variability of the m,’s in its estimation of c. 

The author finds the distribution of R based on four chawrvations when the 
m,’s are not all equal. He then plots the power curves of the range chart for three 
special cases, approximating a two fold integral by using “Circular Probability 
Paper” [1 ]. 

He considers the following three cases: 


Case 1. m; = a + (i — 1)m where a and 7 are constants. 
Case 2. m; = mj = m + 6 = m, + 6 (i, j, k, or 1 represents one of four integers 
1, 2, 3, 4, with no two equal). 
Case 3. mM; = Mm; = m = m, + 6. 
Table 1 gives the probability that a point plotted on a range chart falls outside 
a specified control limit for the three different cases. Table 2 gives the probability 
that x; and x; take on the smallest and largest value for Case 1 where 7 = .05. 
Some typographical errors are listed below: 


1. Fig. 2 Case 1 should be Fig. 3 Case 2. 
2. Fig. 3 Case 2 should be Fig. 2 Case 1. 
3. In Table 2, the headings Max. and Min. should be interchanged. 


SEYMOUR GEISSER 
National of Standards 
25, D. C. 


. F. C. Lone & C. W. Topp, “Circular probability paper,” Industrial Quality Control, v. 9, 
1952, p. 10-16. 


25[K ].—C. T. Fan, “‘Note on the construction of an item analysis table for the 
high-low-27-per-cent method,” Psychometrika, v. 19, 1954, p. 231-237. 


The table whose construction is described in this article provides a means of 
translating the observed proportions of success in the high and low 27 percent 
groups (denoted by pz and pz) into measures of item difficulty (denoted by p 
and A) and of item discrimination (denoted by r). Values are tabled for the diffi- 
culty index, p, and the discrimination index, r, as functions of pz and pz. These 
results are based on Kart PEarsON’s tables of the normal bivariate surface. 
A second item difficulty index, A, can be determined from # and is expressed in 
terms of a normal curve deviate with mean 13 and standard deviation 4. The 
final complete table is not given in this paper but has been published by the 
Educational Testing Service, Princeton, N. J. However, a graphical version of 
the final results is given as figure 2 in this paper. Preliminary computations in the 
table construction are presented in table 3, which furnishes py and pz values to 
4D for .5000 < p < .9713 and r = .05(.05)1.00. Figure 2 furnishes a simplified 
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version of the final chart for use in graphical estimation of p and r for given values 
of px and pz. This graph contains a square grid for p, r = .00(.10)1.00 and curves 
of corresponding pq and py, values for pxpax = .10(.10).90, .99 and p, = .01, 
-10(.10).90. 


Joun E. WatsH 
Military tions Research Division 
Lockheed Aircraft Corporation 
Burbank, California 


26[K ].—W. T. Suarp, J. M. Kennepy, B. J. Sears, & M. G. Hoyze, Tables of 
Coefficients for Angular Distribution Analysis, Atomic Energy of Canada, Ltd., 
Chalk River, Ontario, A.E.C.L. Report No. 97, 1954, xxxix + 38 p., $2.00 
The functions of angular momentum quantum numbers introduced by Racah 

have related closely the studies of the various correlations in angle of particles 
and radiation absorbed or emitted by nuclei. In the Chalk River tables the authors 
have compiled and calculated values of the Racah W-function, the associated 
coefficients Z and Z;, and the 9-7 symbol, X, for ranges of parameters which occur 
in low energy nuclear reactions. 

The coefficients W(ijl'j’; sk) and Z(Jjl'j’; sk) are in the main taken from the 
tables of BIEDENHARN [1] and Ost et al [2]. The range of s, the channel spin 
parameter, is extended to 7/2; a few values for s = 4 not given by Obi are newly 
computed. Except in a few cases, j’ = j; the range of 7 is 0(4)5. / and I’ are 
restricted to integers less than 4, whether or not the selection rules give non- 
zero coefficients for higher integers. A table of W(jj:jj:; Lk) is provided for 
= and L = 1, 2. 

The coefficients Z,(LjL’j’; Ik) are for j’ = 7 closely related to the F,(LI}j) 
and G,(LL’Ij) of BrEDENHARN and Rose [3]. The present tabulation, for L, 
L’ = 1, 2, 3, 7’ = j = 0($)6, I = 0($)6, extends previous tables [4] to the 
parameter values J = 5, 6 and j = 6. 

In all cases & is restricted to the even integers. 

abc 

New tables of the coefficient X|def| are given for a, d = 1,2,e = b, f = ¢, 

ghk 

b, c = 1($)5 and g, h, k = 0(2)8. A second specialization is made to k = 1, 

a = d = 1(3)3, g = h = 2, 4, with the other parameters in the range 1(4)3. 
Finally, auxiliary tables are provided for the triangle coefficient A(a, 6, c), and 

for the CLEBSCH-GORDON coefficients 1, /’ = 0(1)6, and (LL’ — 11|k0), 

L, L’ = 1, 2, 3. 

The square of a coefficient is in each case a rational fraction. The quantities 
tabulated are the (positive or negative) powers of the prime factors of each 
square, and the sign of the coefficient. 

A useful introductory section lists many of the important properties of the 
coefficients and discusses their application in the theory of nuclear reactions, 
heavy particle-gamma ray correlation, gamma-gamma correlation, transforma- 
tions between angular momentum coupling schemes, triple correlation, and 
transition probabilities for electromagnetic radiation. 
As P. V. C. Houca 
Ann Arbor, 
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C. BreDENHARN, Oak Ridge National Laboratory, Reports No. 1098 and 1501. 
Ost et al., Annals of the Tokyo Astronomical Gteervelory III, 1953, p. 89. 
C. Brepennarn & M. E. Rose, Rev. Mod. Phys., . _ 1953, p. 729. 

T. SHarp, Gove, & Taut, Chalk River Reports. PD. no. 254, unpublished. 


y 
3. 
4. W. 
27[K].—J. M. Kennepy, B. J. Sears, & W. T. SHarp, Tables of X Coefficients, 
Atomic Energy of Canada, Ltd., Chalk River, Ontario, 1954, A.E.C.L. 
Report No. 106, ix + 16 p., $1.00. 
This work extends the authors’ table of X coefficients in ‘Tables of Coeffi- 
cients for Angular Distribution Analysis” [1]. The following tables are provided: 


abc 
1. X*\abc| for a = 1, 2, b, c = 1(1)5. Here and throughout g, h, k run over all 
ghk 
even integers compatible with the selection rules for X. 


1bc 
2. X*\2bc| for b, c = 1(1)5. 
ghk 


abe 
3. X*labc| for a = 1, 2, b, c = $(1)$; and for a = 3, 4, b = 4(1)$, c = 3(1)§. 
ghk 


1bc 
4. X*\2bc| for b, c = $(1)§. 
ghk 


akc 
5. X*\abc| for b = 3, $, a = 1(1)4, c = 34(1)§. 
ghk é 


P. V. C. 
University of Michigan 
Ann Arbor, Michigan 


1. W. T. Saarp, J. M. KENNEDY, B. J. SEars, & M. G. Hoye, A.E.C.L. Report No. 97, 1954, 
Chalk River, Ont. See preceding review. 


28[L].—N. W. McLacuaian, Bessel Functions for Engineers, Oxford University 
Press, New York, 1955. Second edition. xii + 239 pages, 24 cm. Price, $5.60. 


This second edition is considerably larger than the first edition, which ap- 
peared about twenty years ago. It still is a satisfactory text on the use of Bessel 
functions in classical physics and engineering problems. A little more has been 
added on the mathematical properties of the functions and a number of different 
applications have been added. The tone is strictly practical throughout, funda- 
mental equations are quoted, discussed and applied, but not always proved. The 
style is condensed but readable. Short five-place tables, included at the back of 
the book, are for J,(x) for n = 0, 1, x = 0(0.1)16; for m = 2, 3, 4, x = 0(0.1)5; 
for zeros of these J’s; for Y,(x) and H,(x), = 0,1, x = 0(0.1)16; for Io, I:, Ko, 
K, for x = 0(0.1)10; and for similar ranges of Ber and Bei, Ker and Kei, and 
their derivatives. 


Puitrp M. Morse 
Massachusetts Institute of Technology 
Cambridge, Massachusetts 
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29[L ].—Monrio Onoe, “Formulae and tables, the modified quotients of cylinder 
functions,” Report of the Institute of Industrial Science, University of Tokyo, 
v. 4, 1955, No. 5 (Serial No. 32), 22 p. 


C,(z) being any cylinder function [1], the author sets 
€,(z) = €,(s) = 


and similarly for Bessel functions of the first, second, and third kinds. He gives 
a collection of formulas regarding these functions, and’ also numerical tables to 
4 to 7S, each accompanied by graphs. 

Table 1. $,.(x) for m = 1(1)5, x = 0(.1)5 and the first zero and the first pole 
of $1(x) with corresponding values of %,. 

Table 2. $,(¢y) for nm = 1(1)5, x = 0(.2)6. 

Table 3. 9),(x) for m = 1(1)5, x = 0(.1)5 with the first two zeros of 9), ail 
the first zeros of Jeo, Js, Da. 

Table 4, 5, (iy), m = 1(1)5, x = 0(.1)5. 

Table 5. $,(”), 9.(m) for m = 1(1)50. 

Appendix. Some coefficients in power series expansions of {{, and asymptotic 
expansions of 5,, m = 1(1)5. 

A, Erpéy1 


California Institute of Technology 
Pasadena, California 


1. G. N. Watson, Bessel Functions, Cambridge, 1922, p. 82. 


30[L ].—H. E. Sauzer, ‘Complex zeros of the error function,” J. Franklin Insti- 
tute, 260, 1955, p. 209-211. 


Let f(z) = So edu. The zeros z, of f(z), for which 44 < argz, < $x, and 
0 <|2:| <|z2|<--- are tabulated as follows: m = 1 to 12D, m = 2(1)10 to 9D. 
This table extends and corrects one of T. LarsLe (Z. ang. Math. u. Phys., 1951, 
p. 484-486). 


31[L].—Haro_p OsTERBERG & Gorpon L. WaLKER, Table of So [J:(x)/x 
Research Center, American Optical Company, Southbridge, Massachusetts, 
Communication No. 1, September 1, 1955. 7 multilithed pages with cover. 
28 cm., deposited in the UMT FILe. 


This table lists the function f(z) = Jo? xJi(x)dx for z = .01(.01)3.85 and 
z = 4(1)25 plus about seven special values of the argument. This is the function 
g(z) of “A guide to tables of Bessel Functions,” by R. C. Arcnipatp, MTAC, 
v. 1, 1944, p. 247. The numbers shown are usually. 7S, and accuracy is claimed 
through 5S or 6S. 

For the early range, the calculations were carried out on an International 
Business Machines Card Programmed Calculator, and for the later range they 
were carried out on a desk machine, using tables of Bessel Functions already 
published. 

A similar, less extensive table has been published by R. Gans, “Mikro- 
skopische probleme,” Annalen der Physik, s. 4, v. 78, 1925, p. 1-34. The authors 
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note that their function is equal to /o* Jo(x)dx — Ji(z). The function J" Jo(x)dx 
has been tabulated to 10D by A. N. Lowan and M. ABramow!7z in “Tables of 
integrals of Jo(x) and Yo(x),” J. Math. and Phys., v. 22, 1943, p. 2-12, and this 
table has been reprinted as AMS 37 by the National Bureau of Standards; its 
range is <=0(.01)10. 


32[L].—B. Zonpek, ‘“‘The values of and ['(3) and their logarithms accurate 
to 28 decimals,”” MTAC, v. 9, 1955, p. 24-25. 


TABLE ERRATA 
Reviews in this issue mention errata in the following works: 


THE RAND CorporaTION, One Million Digits and 100,000 Normal Deviates, 
Review 11, p. 39-43. 

BENJAMIN EpsTEIN, “Truncated life tests in the exponential test,” Ann. Math. 
Stat., v. 25, 1954, p. 555-564, Review 15, p. 44-45. 

R. A. BrapLey & M. E. Terry, “Rank analysis of incomplete block designs. I. 
The method of paired comparisons,” Biometrika, v. 39, 1952, p. 324-345, 
[MTAC, v. 8, 1954, p. 17], Review 22, p. 49. 

SHozo SuimapA, ‘“‘Power of R-charts,” Reports of Statistical A pplication Research, 
Union of Japanese Scientists and Engineers, v. 3, 1954, p. 70-74, Review 24, 
p. 50. 

T. LarBie, “Héhenkarte des Fehler-integrals,” Z. ang. Math. u. Phys., 1951, 
p. 484-486, Review 30, p. 53. 


247.—GIUSEPPE PALAMA & L. PoLetti, ‘‘Tavola dei numeri primi dell’intervallo 
12 012 000-12 072 060,’’ Unione Matematica Italiana, Bollettino, s. 3, v. 8, 
1953, p. 52-58. (MTAC, v. 7, 1953, p. 173, Review 1101[F].) 

The following errata have been found. 


Entry Division Probably intended prime 
12 019 307 277 ~ —_ 
12 020 023 1901 = 
12 023 381 31 12 023 383 
12 028 813 131 12 028 817 
12 045 149 457 — 
12 047 023 107 — 
12 071 881 2081 “= 


In addition the following primes should be added to the list. 


12 047 309 
12 069 919 
N. G. W. H. BEEGER 


Holland 


Editor's note: Primality of each number listed above as prime has been verified on the SWAC 
computer by J. L. SELFRIDGE. 


CORRIGENDA 


CORRIGENDA 


Review 56[C, D, E, K, L, S].—Crci Hastincs, jr., JEANNE T. Haywarp, & 
James P. Wonca, jr., “Approximations for Digital Computers,” MTAC, v. 9, 
p. 122, line 6, for k; read K, in 


du 
Also replace the sentence, ‘Many of these approximations have appeared 


previously in MTAC, v. 7 and v. 8” by “Many such approximations have ap- 
peared in MTAC, v. 7 and v. 8.” 
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